“Suspensões de Veículos – Um Estudo de Caso” MILTON KOITI AKIYAMA 259 UNESP Faculdade de Engenharia do Campus de Guaratinguetá Guaratinguetá 2005 MILTON KOITI AKIYAMA SUSPENSÕES DE VEÍCULOS - UM ESTUDO DE CASO Dissertação apresentada à Faculdade de Engenharia do Campus de Guaratinguetá, Universidade Estadual Paulista, para a obtenção do título de Mestre em Engenharia Mecânica na área de Projetos e Materiais. Orientador: Prof. Dr. José Elias Tomazini Guaratinguetá 2005 A315s Akiyama, Milton Koiti Suspensão de veículos – um estudo de caso / Milton Koiti Akiyama.- Guaratinguetá : [s.n.], 2005 66f.:il. Bibliografia: f. 60-61 Dissertação (mestrado) – Universidade Estadual Paulista, Faculdade de Engenharia de Guaratinguetá, 2005 Orientador: Prof. Dr. José Elias Tomazini 1.Veículo 2.Suspensão passiva 3.Rugosidade de pista I. Título CDU629.015 DADOS CURRICULARES MILTON KOITI AKIYAMA NASCIMENTO 21.03.1946 – SÃO PAULO / SP FILIAÇÃO Shigeo Akiyama Fussae Akiyama 1966/1971 Curso de Graduação Engenharia Mecânica – Escola de Engenharia de Taubaté 2003/04 Curso de Pós-Graduação em Engenharia Mecânica, nível de Mestrado, na Faculdade de Engenharia do Campus de Guaratinguetá da UNESP. AGRADECIMENTOS Em primeiro lugar agradeço a Deus, fonte da vida e da graça. Agradeço também a todas as outras entidades que me deram alento, força, conselhos, à minha família e meus amigos. Ao meu orientador, Prof. Dr. José Elias Tomazini, que sempre me orientou e incentivou, teve paciência, mesmo nos meus momentos de desânimo. Sem a sua orientação, auxílio, disposição e dedicação, o estudo aqui apresentado teria sido praticamente impossível. Ao meu pai, Shigeo, que Deus o tenha, que sempre me incentivou, e me mostrou com o seu próprio exemplo de tenacidade, o caminho a seguir; À minha mãe, Fussae, que sempre me apoiou com palavras doces os esforços. À minha esposa, Midori, que nunca deixou que eu me abatesse, dando- me apoio espiritual e palavras de estímulo. E também gostaria de agradecer ao pessoal da Secretaria de Pós- Graduação, em especial a Regina Célia Galvão Faria Alves e Elisa Mara de Carvalho Nunes, que sempre demonstraram boa vontade e simpatia, esclarecendo e orientando. AKIYAMA, M.K. - Suspensões de Veículos – Um Estudo de Caso. 2005. Dissertação (Mestrado em Engenharia Mecânica) - Faculdade de Engenharia do Campus de Guaratinguetá, Universidade Estadual Paulista, Guaratinguetá, 2005. RESUMO Neste trabalho é feito um estudo de comportamento dinâmico de um veículo pequeno com suspensão passiva trafegando em pista com ondulações regulares, e análise da estabilidade. O veículo é modelado considerando-o como um sistema multicorpo rígido, com sete graus de liberdade, utilizando-se o método de Kane. As equações do movimento são definidas sobre uma listagem de programa preparada segundo o programa AUTOLEV, o qual emite um programa em linguagem FORTRAN para as simulações. Foi necessário obter-se diretamente do veículo dados característicos como: centros de gravidade, coeficientes de rigidez dos pneus e das molas, e coeficientes de amortecimento para a realização das simulações. Além disso, foram definidos dois tipos de rugosidade, quatro velocidades, duas seqüências de rugosidade, sendo uma com cinco ondas e outra com vinte ondas, bem como a consideração do veículo com um ocupante e cinco ocupantes. Com os resultados do processamento dos dados, foram traçados gráficos e realizadas análises dos efeitos das rugosidades sobre a dinâmica do veículo, conclusões e apresentações de sugestões. PALAVRAS-CHAVE: Veículo, suspensão passiva, estabilidade, rugosidade de pista AKIYAMA, M.K. Suspensions of Vehicles – A Case Study 2005 Dissertação (Mestrado em Engenharia Mecânica) - Faculdade de Engenharia do Campus de Guaratinguetá, Universidade Estadual Paulista, Guaratinguetá, 2005. ABSTRACT In this paper we study the dynamic behaviour of a small vehicle provided with passive suspension, riding on regularly rugged roads, and stability analysis. The vehicle is modelled considering as a rigid multibody system, with seven degrees of freedom, applying the Kane´s method. Motion equations are defined on a program list according to the AUTOLEV software, which gives a FORTRAN language program for simulations. It was necessary to obtain directly from the vehicle, characteristic datas as centres of gravity, tires and springs stiffnesses, damper characteristics. In addition, it was defined two road roughnesses, four riding speeds, two sequences of roughness, one with five and other with twenty waves, as well as the consideration of vehicle with one and five occupants. With the results of processing, graphics were plotted, and analysis realized about roughness effects on vehicle dynamics and presentations of suggestions. KEYWORDS: Vehicle, passive suspension, stability, road roughness LISTA DE FIGURAS Figura 1 - Rugosidade de pista 20 Figura 2 - Representação de um sistema de suspensão simples 22 Figura 3 - Modelo do veículo considerado, com sete graus de liberdade 28 Figura 4 - Seqüência de passagem do veículo pela ondulação 30 Figura 5 - Determinação da posição do CG do veículo vazio 31 Figura 6 - Veículo erguido para medição da carga sobre o eixo dianteiro 32 Figura 7 - Numeração dos pontos de referência no veículo – Vista em Elevação 34 Figura 8 - Numeração dos pontos de referência no veículo – Vista em Planta 34 Figura 9 - Carregamento do veículo com passageiros – Vista em Elevação 34 Figura 10 - Carregamento do veículo com passageiros – Vista em Planta 35 Figura 11 - Determinação das reações de apoio dianteiro e traseiro com o veículo carregado 35 Figura 12 - Esquema para determinação do coeficiente de rigidez do pneu 36 Figura 13 - Esquema para determinação do coeficiente de rigidez da mola de suspensão principal 39 Figura 14 - Montagem esquemática da suspensão da roda dianteira 40 Figura 15 - Esquema da suspensão traseira com o veículo vazio e carregado 41 Figura 16 - Gráfico DYA10- 501 A 504 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 46 Figura 17 - Gráfico DYA10- 505 A 508 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos, veículo com cinco pessoas. 47 Figura 18 - Gráfico DYA10- 601 A 604 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 48 Figura 19 - Gráfico DYA11- 501 A 504 – Oscilação vertical das rodas traseiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 48 Figura 20 - Gráfico DYA11- 601 A 604 – Oscilação vertical das rodas traseiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 49 Figura 21 - Gráfico FA1- 501 A 504 – Força na conexão superior do amortecedor dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos 50 Figura 22 - Gráfico FA1- 601 A 604 – Força na conexão superior do amortecedor traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 51 Figura 23 - Gráfico FA3- 501 A 504 – Força na conexão superior do amortecedor dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 52 Figura 24 - Gráfico FA3- 601 A 604 – Força na conexão superior do amortecedor traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 52 Figura 25 - Gráfico TORQ-501 A 504 – Torque no braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 53 Figura 26 - Gráfico TORQ- 601 A 604 – Torque no braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 54 Figura 27 Gráfico Q4- 501 A 504 – Rotação do braço de suspensão dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos 55 Figura 28 - Gráfico Q4-601 A 604 – Rotação do braço de suspensão dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 55 Figura 29 - Gráfico Q6-501 A 504 – Rotação do braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 56 Figura 30 - Gráfico Q6-601 A 604 – Rotação do braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 56 LISTA DE TABELAS Tabela 1 Deflexões das suspensões e rodas do Clio 37 Tabela 2 Tempo de percurso para trecho com 5 ondulações 45 Tabela 3 Tempo de percurso para trecho com 20 ondulações 46 ‘ LISTA DE ABREVIATURAS E SIGLAS MIRA British Motor Industry Research Association NBR Norma Brasileira Registrada SI Sistema Internacional LISTA DE SÍMBOLOS CD Coeficiente de amortecimento do amortecedor dianteiro N.s/m CT Coeficiente de amortecimento do amortecedor traseiro N.s/m df1 a df4 Deflexão verificada na roda do veículo, quando o veículo é carregado m df5 a df8 Deflexão verificada na marca sobre a carroçaria do veículo, quando o mesmo é carregado m dff Deflexão do pneu dianteiro adotado com o veículo carregado m dft Deflexão do pneu traseiro adotado com o veículo carregado m fv1 a fv4 Distância vertical medida do centro das rodas 1 a 4 até o solo, veículo vazio m fv5 a fv8 Distância vertical medida das marcas 5 a 8 no veículo até o solo, veículo vazio m fk1 a fk4 Distância vertical medida do centro das rodas 1 a 4 até o solo, veículo carregado m fk5 a fk8 Distância vertical medida das marcas 5 a 8 no veículo até o solo, veículo carregado m Hgv Altura do Centro de gravidade do veículo vazio m KFapp Coeficiente de rigidez aparente da mola de suspensão dianteira, no plano da roda N/m KF Coeficiente de rigidez efetivo da mola principal dianteira do veículo N/m KT Coeficiente de rigidez à torção da mola principal traseira do veículo N.m/rad KPf Coeficiente de rigidez do pneu dianteiro N/m KPt Coeficiente de rigidez do pneu traseiro N/m L Distância entre eixos do veículo m Lfv Distância do CG do veículo vazio ao eixo frontal m Lfk Distância do CG do veículo carregado ao eixo frontal m Ltv Distância do CG do veículo vazio ao eixo traseiro m Ltk Distância do CG do veículo carregado ao eixo traseiro m Lfa,b Distância das cargas Qa e Qb em relação ao eixo frontal do veículo m Lfc,d,e Distância das cargas Qc, Qd e Qe em relação ao eixo frontal do veículo m n Altura de levantamento das rodas traseiras para determinação da altura do CG do veículo m Qa,b Cargas Qa e Qb, no mesmo alinhamento, na figura N Qc,d,e Cargas Qc, Qd e Qe no mesmo alinhamento, na figura N Rf Raio de rolamento da roda dianteira m Sf Bitola frontal do veículo m St Bitola traseira do veículo m Wv Peso do veículo vazio N Wfv Carga das rodas dianteiras ao solo, com o veículo vazio N Wfv1 Carga de uma roda dianteira, com o veículo vazio N Wtv Carga das rodas traseiras ao solo, com o veículo vazio N Wtv1 Carga de uma roda traseira, com o veículo vazio N W´f Carga do eixo frontal sobre o solo, com o eixo traseiro levantado N Wk Peso total do veículo, carregado N Wfk Carga de apoio das rodas frontais ao solo, com o veículo carregado N Wfk1 Carga de uma roda dianteira, com o veículo carregado N Wtk Carga de apoio das rodas traseiras ao solo, com o veículo carregado N Wtk1 Carga de uma roda traseira, com o veículo carregado N ∆Wt Variação de carga sobre uma roda traseira, para veículo vazio e carregado N ∆θ Rotação do braço de suspensão grau SUMÁRIO LISTA DE FIGURAS LISTA DE TABELAS LISTA DE ABREVIATURAS E SIGLAS LISTA DE SÍMBOLOS 1 INTRODUÇÃO 18 2 CONCEITOS PRELIMINARES 20 2.1 RUGOSIDADE DE PISTA 20 2.2 FUNÇÃO DOS ELEMENTOS DE SUSPENSÃO 22 2.2.1 Molas 23 2.2.2 Amortecedores 24 2.3 MODELAGEM DO VEÍCULO PELO MÉTODO DE KANE – CONCEITOS BÁSICOS 25 2.4 MODELAGEM DO VEÍCULO 27 2.5 USANDO O PROGRAMA AUTOLEV 29 3 DETERMINAÇÃO DE PARÂMETROS PARA MODELAGEM DO VEÍCULO 31 3.1 DETERMINAÇÃO DO CENTRO DE GRAVIDADE DO VEÍCULO VAZIO 31 3.2 DETERMINAÇÃO DOS COEFICIENTES DE RIGIDEZ 33 3.2.1 Determinação do Coeficiente de Rigidez dos Pneus 36 3.2.2 Determinação do Coeficiente de Rigidez das Molas de Suspensão 38 3.2.2.1 Coeficiente de Rigidez das Molas de Suspensão Frontal 39 3.2.2.2 Coeficiente de Rigidez das Molas de Suspensão Traseira 40 3.2.2.3 Dados dos Amortecedores 42 4 SIMULAÇÕES 43 4.1 PROCESSAMENTO NO PROGRAMA FORTRAN 43 4.2 PLOTAGEM DAS CURVAS 43 4.3 CONDIÇÕES DE TRÁFEGO CONSIDERADAS NAS SIMULAÇÕES 43 4.4 RESULTADOS DE SIMULAÇÕES 44 4.5 ANÁLISE DOS RESULTADOS 45 5 CONCLUSÕES 58 REFERÊNCIAS 61 APÊNDICE 1 – Listagem do programa AUTOLEV para a versão Clio 6.3 APÊNDICE 2 – Listagem do programa FORTRAN para a versão Clio 6.3 APÊNDICE 3 – Dados de entrada das simulações efetuadas APÊNDICE 4 – Gráficos dos ensaios 1. INTRODUÇÃO Embora o sistema de suspensão das rodas de um veículo pareça simples em concepção quando comparado com o motor ou a transmissão, a sua função e o seu efeito sobre o comportamento do veículo são muito complexos. A grande maioria dos estudos atuais está dirigida aos veículos providos de suspensões ativas e semi-ativas, cujo desenvolvimento é resultado das aplicações em veículos de competição tipo Fórmula 1 e Formula Indy, entretanto o custo final para os usuários comuns ainda é muito alto. O presente trabalho tem por objetivo estudar o comportamento de um pequeno veículo urbano com suspensão passiva trafegando em pista com rugosidades a intervalos regulares e a velocidades diversas, para análise de estabilidade do veículo e de conforto dos ocupantes, e propor soluções. A principal preocupação nesses estudos está em verificar se em alguma condição de tráfego pode ocorrer ressonância na suspensão gerando instabilidade, perda de contato com o solo, esforço acima dos valores admissíveis caso o movimento da suspensão ultrapasse os cursos dos amortecedores ou das molas, e outras conseqüências. Com relação a estudo sobre suspensões passivas, pode-se citar o artigo de SHARP, R.S. e HASSAN, S.A. (Sharp, 1986), onde um quarto de veículo é analisado com a sua suspensão submetida a perturbações aleatórias, sendo definidos alguns valores diferentes de rigidez e de amortecimento, e posteriormente os resultados são discutidos e apresentadas algumas sugestões. Outro artigo que aborda a questão das suspensões de veículos é de SHARP, R.S. e CROLLA, D.A. (Sharp, 1987), onde as suspensões de veículos, ativas, semi-ativas e passivas são analisadas em condições de tráfego sobre pistas com rugosidade aleatória. No caso de veículos com suspensão passiva, a análise é realizada em um quarto de veículo e o artigo menciona o tradicional conflito entre o conforto e o curso de trabalho da suspensão: normalmente, para uma 18 mesma rugosidade de pista e velocidade de translação do veículo, o conforto pode ser aumentado diminuindo a rigidez e o amortecimento; enquanto isso, o espaço de trabalho da suspensão aumenta muito. Para o presente estudo são utilizados os dados obtidos a partir de um Clio 1.6, ano de fabricação 1998 de propriedade da família, fabricado na Argentina. Para a realização de ensaios em computador, o veículo utilizado como referência foi considerado como um modelo multicorpo rígido, adotando-se o método de Kane para modelagem, conforme comentado em linhas gerais no Capítulo 2, onde são apresentados conceitos básicos sobre os modelos multicorpos rígidos, as vantagens da aplicação do método de Kane, e os programas que facilitam a aplicação do método. O modelo adotado apresenta sete graus de liberdade, escolhidos para se obter o maior número de informações sobre o comportamento do veículo. Para adequar ao processamento, durante a modelagem foram considerada s algumas alterações nas características do veículo, procurando-se desviar o mínimo possível em relação à configuração original do veículo. Os tipos de rugosidade de pista mencionados nas referências são relacionados no capítulo 2, sendo que, para fins de estudo apenas dois tipos foram adotados. No item 2.4, é feita uma descrição da modelagem do veículo, incluindo como foi decidido utilizar-se sete graus de liberdade. No item 2.5, é feita uma descrição do uso do programa AUTOLEV que efetua os cálculos necessários, baseado no método de Kane e cria um programa na linguagem FORTRAN para as simulações. Os resultados do processamento são apresentados utilizando-se um programa para plotagem de gráficos. No capítulo 3, estão descritos os procedimentos experimentais efetuados para a obtenção de dados do veículo para inserção nas simulações, como a localização do Centro de Gravidade do veículo quando vazio e quando cheio, as características de rigidez dos pneus, das molas de suspensão dianteira e traseira, e os dados dos amortecedores dianteiro e traseiro. Nesse processo, foi necessário 19 realizar algumas aproximações, pois as medições apresentaram algumas diferenças, provavelmente devido a alguns componentes gastos no veículo. No capítulo 4, definem-se as variáveis utilizadas nas simulações, como velocidades, tipos de rugosidade selecionados, quantidades de ondas percorridas e condições de carga do veículo. Os principais gráficos obtidos nas simulações são analisados no item 4.5, procurando-se interpretar as causas dos comportamentos demonstrados nos mesmos. No capitulo 5, são apresentadas as conclusões finais e sugestões de estudos com o objetivo de melhorar as respostas das suspensões do veículo. Com a finalidade de contribuir para o desenvolvimento de novos trabalhos de pesquisa, apresentam-se algumas sugestões para futuros trabalhos. 20 2 CONCEITOS PRELIMINARES 2.1 RUGOSIDADE DE PISTA Nas referências são descritos diversos tipos de irregularidades encontradas nos pavimentos pelo mundo, desde as mais castigadas pela falta de manutenção, mas também pela execução sem a preparação adequada do contra-piso, e as pistas em terra ou areia que apresentam ondulações produzidas pelo tráfego intenso de veículos pesados. Figura .1 – Rugosidade de pista Segundo Bastow (1993), de um modo geral as estradas que apresentam amplitudes inferiores a 0,005 m (Fig. 1) são consideradas muito boas, e são de média qualidade as que expressam amplitudes máximas de 0,013 m. As estradas muito boas com amplitude de 0,005 m têm passo entre 4 e 5 m, com taxa acumulada de 1,18 m/km, representando 236 a 250 calombos por quilometro. Por outro lado, as ondulações superficiais com amplitudes de 0.019 a 0,025 m sempre causam desconforto nas velocidades normais, da ordem de 80 a 100 km/h. Nas literaturas analisadas são mencionadas as seguintes rugosidades: • Pavimento belga – representa um tipo de pista encontrada em regiões da Bélgica, França, Holanda e Alemanha durante a Segunda Guerra, com blocos de pedra irregulares, sendo que há uma pista de testes construída pela MIRA (British Motor Industry Research Association) representando uma seção típica da Bélgica com pedras arredondados de amplitude de 0,025 m e passos aleatórios entre 0,15 e 0,23 m, utilizada para testes de resistência do veículo. É consenso geral que se 21 um protótipo sobrevive depois de submetido a esse tipo de pavimento durante 1600 km, é improvável que desenvolva defeitos estruturais durante o período de vida normal (BASTOW, 1993). Esse tipo de pista não é considerado no presente trabalho. • Pavimentos de grande deterioração - encontrados em paíises com manutenção mínima e condições climáticas que produzem deterioração severa. Os fabricantes de veículos que pretendem desenvolver e manter mercado para seus produtos nessas regiões devem projetar e desenvolver suspensões capazes de resistir ao uso continuo nessas condições. A pista de testes construída pela MIRA com representação desse tipo tem rugosidade com amplitude de 0,025 m e passos de 0,75 m + 0,05 m. • Pista com longas ondulações – Neste caso as superfícies apresentam ondas senoidais com amplitude de 0,1 m e 12 m de passo. Esse tipo de pista tende a excitar a massa suspensa com 1,3 a 1,5 Hz, a velocidades entre 55 e 65 km/h. Há dois motivos para se testar veículos nessas condições: primeiro, como representa uma pista nova, os fabricantes querem assegurar que seus veículos podem trafegar a velocidades razoáveis em tais superfícies e segundo, permite testes de durabilidade em veículos de série aplicando grande amplitude na suspensão. Esse tipo de ondulação também não é considerado neste trabalho. • Estrada de terra batida– amplitude de 0,015 a 0,02 m e passos entre 0,50 e 0,60 m. • Estradas precárias, clima severo – amplitudes de 0,050 m e passos de 0,75 m; • Estradas precárias em terra batida– amplitude de 0,015 m e passos de 0,15 a 0,20 m; • Estradas precárias em areia fofa– amplitude de 0,025 m e passos de 0,25 a 0,30 m 22 2.2 FUNÇÃO DOS ELEMENTOS DE SUSPENSÃO A função primária do sistema de suspensão de um veículo é isolar a estrutura (e os seus ocupantes) de choques e vibrações geradas pelas irregularidades da superfície da estrada, de modo a conciliar a sensibilidade do ser humano, enquanto mantém a estabilidade, o controle direcional e todas as necessidades de manobra do veículo em seu comportamento dinâmico. A sensibilidade humana às vibrações é bastante complexa, sendo estudada com detalhes por Donald Bastow (BASTOW, 1993). De um modo geral, freqüências de vibrações verticais entre 1,5 e 2,3 Hz são consideradas confortáveis, bem como são aceitáveis as oscilações longitudinais ou laterais abaixo de 1,5 Hz. Tonturas e enjôos são provocados quando o ouvido interno é submetido a baixas freqüências, entre 0,5 e 0,75 Hz, e sérios desconfortos são sentidos no organismo a freqüências de 5 a 7 Hz. Entretanto, a função básica de um sistema de suspensão em um veículo não é proporcionar conforto aos ocupantes quanto à vibração, embora seja desejável, mas sim à manutenção do contato entre as rodas e a superfície da estrada, pois o controle direcional e a estabilidade do veículo dependem disso. Figura 2 - Representação de um sistema de suspensão simples 23 Um sistema de suspensão pode ser representado de modo mais simples considerando um quarto de veículo, como mostrado na figura.2. Uma massa suspensa representada pelo corpo do veículo e os seus ocupantes, apoiado por intermédio de uma mola sobre a massa não suspensa representada pelo conjunto da roda, freio, os mecanismos de articulação da suspensão, os braços de rigidez e outros, que por sua vez apoiam-se ao solo com a rigidez do pneu atuando como uma mola. Quando momentaneamente perturbada, a massa suspensa oscila verticalmente com a sua freqüência natural, e continua oscilando devido ä ação das molas da suspensão e do pneu combinadas. Para eliminar rapidamente as oscilações é montado um absorvedor de choques combinado com a mola de suspensão principal. O absorvedor de choques normalmente é um amortecedor hidráulico, que reage com forças de resistência contrárias aos movimentos oscilatórios e proporcionais ao quadrado da velocidade da massa suspensa em relação à massa não suspensa. O pneu também apresenta uma pequena parcela de amortecimento, e não foi considerado no presente estudo. Os atritos internos no sistema de suspensão também contribuem para o amortecimento das oscilações, sendo constantes e independentes da velocidade e da amplitude. Esses atritos também foram desconsiderados no presente estudo. Para isolar o veículo das irregularidades da pista, a suspensão deve agir como um filtro. Entretanto, deve também sustentar o peso estático do veículo, criando assim uma deflexão estática na mola, que deve ser levada em conta no projeto. Deve também levar em conta as variações de carregamento do veículo, sem considerar outra forma de sustentação. (WILLIAMS R.A., 1997) 2.2.1 Molas Um dos elementos principais de um sistema de suspensão é a mola, que é o elemento flexível cuja função é de se defletir quando a roda encontra um ressalto e sofre um rápido impulso para cima (rebound). Se não existisse uma mola entre a roda e a carroçaria, o choque transmitido seria considerável. Com a presença da 24 mola, a intensidade de força transmitida à carroçaria seria a necessária para comprimir a mola suficientemente para que a roda passe pelo obstáculo. Essa força provoca a aceleração do corpo para cima, mas a taxas muito inferiores da que ocorreria se não existisse a mola. Por outro lado, quando a roda cai em uma depressão, a mola que está pré- comprimida atua sobre a massa não suspensa relativamente menor, ou seja, o conjunto da roda e eixo, e empurra-a para baixo rapidamente de modo que a roda encontra o fundo da depressão quase sempre antes que a massa muito maior da carroçaria, devido a sua inércia, tenha tempo para iniciar a descida. Devido ao fato que as variações de força na mola para tais deflexões são relativamente pequenas, a aceleração da carroçaria para baixo, suportada pela mola é correspondentemente menor se comparado com a aceleração produzida pela ação da gravidade na ausência da mola. Após a passagem pelo distúrbio, seja um ressalto ou uma depressão, o movimento subseqüente da carroçaria é a sua vibração livre sobre as molas, sendo pequena a aceleração. Durante uma deflexão, a mola armazena energia potencial da carga aplicada ao se deformar, seja linearmente ou angularmente, conforme a configuração da mesma. A rigidez da mola é medida como a taxa de carga que produz a deflexão da mola. A unidade usual é Newton por metro (N/m) no SI - Sistema Internacional, considerando as molas por deflexão linear, e de Newton- metro por radiano (N.m/rad), para as molas de deflexão angular. 2.2.2 Amortecedores As molas são combinadas com elementos amortecedores, cuja função é de limitar e amortecer rapidamente as oscilações indesejáveis da massa suspensa do veículo e também para evitar que as rodas percam o contato com o solo. Os amortecedores ou absorvedores de choque são necessários para produzir um rápido desaparecimento de vibrações na freqüência natural produzidas no sistema de suspensão, seja por esforços randômicos ou periódicos e que pode introduzir um estado de ressonância. Para isso, o amortecedor reage 25 com uma força oposta ao movimento instantâneo da suspensão e proporcionais à velocidade da massa suspensa em relação à massa não suspensa sendo que, para pequenas amplitudes e baixa velocidade de oscilação, o efeito de amortecimento é desprezível.. 2.3 MODELAGEM DO VEÍCULO PELO MÉTODO DE KANE – CONCEITOS BÁSICOS Para o estudo da dinâmica do veículo, o mesmo foi representado por um modelo rígido multicorpo, através do método de Kane, para determinar os deslocamentos e esforços. O modelo mais simples que pode-se adotar para o estudo da dinâmica de corpos é a partícula, quando o interesse se resume no movimento do corpo como um todo e as rotações ao redor do seu centro de massa são desconsideradas. Uma extensão deste modelo é o corpo rígido, formado por muitas partículas onde as distâncias entre as mesmas permanecem constantes, mas nesse caso, além do interesse pelo movimento global do corpo, deve-se conhecer o movimento de cada partícula do corpo em relação ao centro de massa. Além disso, quando há vários corpos interligados, onde o movimento de um é influenciado pelo outro, tem-se uma extensão dos modelos descritos acima e denomina-se o modelo como multicorpo. (TOMAZINI, 1996) Portanto, define-se um modelo multicorpo definido como um conjunto formado por corpos rígidos e flexíveis interligados através de molas, amortecedores, juntas rotativas, de translação, atuadores, etc., ou seja, interligados de alguma forma (BIANCHI, 1985; SHABANA, 1989; SCHIEHLEN, 1985). A união dos corpos pode formar uma cadeia aberta ou podem estar internamente vinculados formando uma cadeia fechada. O sistema pode também ser externamente vinculado através do movimento prescrito de algum dos corpos. Conforme descrito acima, diversos mecanismos, robôs industriais, estruturas espaciais, suspensões de veículos, motores, etc., podem ser analisados através do modelo multicorpo. Os sistemas multicorpos para 26 aplicações diversas foram estudados em trabalhos desenvolvidos por Kane (1961, 1978, 1980, 1983, 1985), Kreuzer (1984), Körtum (1985), Schiehlen, (1977, 1985), e Shabana (1985, 1986, 1989). Nos trabalhos de Kane as equações do movimento são expressas em termos das velocidades generalizadas e coordenadas generalizadas do sistema. O número de equações diferenciais que descrevem o movimento é igual ao número de graus de liberdade do sistema. Kane introduz os conceitos de velocidades generalizadas, velocidades parciais de partículas ou pontos, velocidades angulares parciais de corpos, forças ativas generalizadas e forças de inércia generalizadas. Nesta formulação os esforços reativos são eliminados graças ao conceito de forças ativas generalizadas. Caso seja necessário, estes esforços podem ser obtidos através do rompimento dos vínculos que os produzem. As equações resultantes da aplicação do modelo multicorpo são bastante complexas. Para auxiliar os cálculos, surgiram nos últimos anos diversos programas em computador que executam todo o trabalho árduo de se obter a mão as equações do movimento, formulando-as equações e alguns até mesmo resolvendo-as. Inicialmente, surgiram os softwares numéricos, onde as equações numéricas são geradas para cada conjunto de parâmetros e/ou para cada passo num procedimento de integração. Esses métodos têem como desvantagem a baixa precisão e um maior tempo de processamento. Mais recentemente surgiram os formalismos simbólicos que eliminaram as desvantagens acima. As equações são geradas uma única vez e os parâmetros podem ser alterados durante a fase de integração das equações. O AUTOLEV, desenvolvido por Levinson (1991) é baseado no método de Kane. O usuário, através da utilização dos comandos do software, vai calculando passo a passo todas as quantidades necessárias para a obtenção das equações do movimento. Após obtê-las, o AUTOLEV cria um programa na linguagem FORTRAN para a integração das equações. 27 No presente trabalho foi utilizado o método de Kane para o estudo de sistemas multicorpos constituídos apenas por corpos rígidos. No livro “Dynamics: Theory and Applications” de autoria de T.R.Kane e D.A.Levinson (Kane, 1985) o método é apresentado detalhadamente. O método é também apresentado no capítulo 2 do trabalho “O modelo multicorpo aplicado a um manipulador modelo rígido e flexível” (Tomazini, 1997). 2.4 MODELAGEM DO VEÍCULO A decisão sobre o número de graus de liberdade no modelo foi tomada após análise das diversas opções. O modelo com apenas um grau de liberdade considerando um quarto de veículo foi descartado desde o início, por não oferecer variações nos estudos. A alternativa seguinte foi a de modelar com meio veículo, que resultaria em quatro graus de liberdade, ou seja, a oscilação linear vertical, a oscilação angular frontal segundo eixo transversal (arfagem), e as oscilações verticais nas suspensões dianteira e traseira. Mas essa opção limitaria o estudo do comportamento dinâmico do veículo, não refletindo as oscilações laterais (rolagem) e a guinada. O modelo adotado é de um veículo inteiro com sete graus de liberdade, como mostrado na Figura 3. A suspensão dianteira representa os amortecedores e molas montados verticalmente, embora sejam ligeiramente inclinados em função da geometria McPherson. Na suspensão traseira, que apresenta concepção com braços de arraste e mola de torção com eixo transversal ao veículo e amortecedores montados em posição bastante inclinada, na modelagem considera-se o amortecedor na posição vertical, de modo a ter o mesmo braço de momento em relação ao eixo de articulação do braço. As massas das rodas incluem a própria roda com o aro, pneumático, cubo, sistema de freio, braços de suspensão. O restante considera-se como massa suspensa. 28 No centro de gravidade S do corpo do veículo atribui-se um sistema de coordenadas, com os eixos s1, s2 e s3, bem como um sistema de coordenadas newtoniano, fixo em um ponto fixo O do espaço com os eixos n1, n2 e n3, sendo que o índice 1 refere-se ao sentido longitudinal do veículo e orientado para frente, o índice 2 para o eixo horizontal transversal com o sentido positivo para a direita, e o índice 3 para o eixo vertical, orientado para baixo. Para os elementos rígidos do multicorpo, considera-se um sistema de coordenadas localizado no centro de gravidade de cada elemento, com as mesmas orientações indicadas acima. Figura 3 - Modelo do veículo considerado, com sete graus de liberdade 29 Para fins de análise do comportamento do veículo, consideram-se os seguintes movimentos: • q1 – oscilação vertical do centro de massa do veículo, segundo o eixo s3; • q2 - rotação do corpo do veículo em torno do eixo longitudinal s1, representando o rolamento; • q3 – rotação do corpo do veículo em torno do eixo transversal s2, representando a arfagem; • q4 – rotação do braço de suspensão dianteiro esquerdo; • q5 – rotação do braço de suspensão dianteiro direito; • q6 – rotação do braço de suspensão traseiro esquerdo; • q7 – rotação do braço de suspensão traseiro direito; 2.5 USANDO O PROGRAMA AUTOLEV A modelagem do veículo no programa AUTOLEV baseia-se no modelo adotado (Figura 3), com a digitação das equações básicas de posição e de movimento dos elementos do multicorpo, as condições iniciais, notações dos dados de entrada e de saída, equações cinemáticas, matrizes de transformação, vetores, velocidades, acelerações, forças e torques, em um processador de texto básico, como um bloco de notas, o qual é renomeado com a extensão “.atl” para processamento pelo programa AUTOLEV. A massa suspensa do veículo considera-se concentrada no centro de gravidade, denominado S; A massa do conjunto da roda frontal, incluindo a roda, o cubo, componentes do freio, e braços de suspensão, indica-se como MD; A massa do conjunto da roda traseira, incluindo a roda, o cubo, componentes do freio, e braços de suspensão, indica-se como MT; Considera-se também que, durante a passagem do veículo pela região rugosa da pista, há a seguinte seqüência a ser considerada, conforme Figura 4: 30 Figura 4 - Seqüência de passagem do veículo pela ondulação A – Inicialmente, o veículo está trafegando em pista plana, portanto sem oscilações; B – As rodas frontais encontram a primeira ondulação, enquanto que as rodas traseiras ainda permanecem na pista plana; C – As rodas traseiras encontram a primeira ondulação, enquanto que as rodas frontais estão percorrendo as ondulações seguintes; D – O veículo inteiro está percorrendo as ondulações seqüenciais; E – As rodas dianteiras saem da última ondulação e entram em pista plana, enquanto que as rodas traseiras ainda percorrem pista ondulada; F – O veículo está percorrendo pista plana como um todo. A listagem do programa “Clio63.atl” encontra-se no Apêndice 1, sendo que “63” significa “versão 6.3”. 31 3 DETERMINAÇÃO DE PARÂMETROS PARA MODELAGEM DO VEÍCULO 3.1 DETERMINAÇÃO DO CENTRO DE GRAVIDADE DO VEÍCULO VAZIO Para a determinação da posição do Centro de Gravidade (CG) do veículo (Figura 5), utilizou-se o método da medição da distribuição do peso aos eixos dianteiro e traseiro, e pela somatória de momentos (TABOREK, 1957). Para a medição das cargas, utilizou-se uma balança de caminhões de uma pedreira, que gentilmente realizou as medições. O operador afirma que a precisão da mesma é de 10 kgf, mas não foi possível confirmar. Inicialmente, colocou-se o veículo inteiro sobre a balança, para determinação do seu peso, com o mínimo de combustível para não influenciar as medições. Em seguida, mediu-se a carga sobre o eixo dianteiro, e posteriormente sobre o eixo traseiro. Os resultados, emitidos em cartões impressos da própria máquina são os seguintes: Wv = 1000 kgf = 9810 N Wfv = 640 kgf = 6278 N Wtv = 360 kgf = 8532 N Figura 5 - Determinação da posição do CG do veículo vazio. 32 Para a determinação da altura do CG, apoiou-se o veículo pelas rodas dianteiras sobre a balança, e as rodas traseiras foram levantadas por meio de um macaco hidráulico (tipo jacaré), até uma altura n em relação ao solo, e mediu-se a carga W´f do eixo dianteiro (Figura 6), obtendo-se o valor: W’f = 680 kgf = 6671 N Figura 6 - Veículo erguido para medição da carga sobre o eixo dianteiro As seguintes expressões foram utilizadas (TABOREK, 1957): Wtv LLfv Wv ⋅ = (1) L Lfv Ltv= + (2) 2 2(W f́ Wfv) L L n Hgv Rf Wv n − ⋅ ⋅ − = + ⋅ (3) Após os cálculos, a posição longitudinal do CG foi determinada com os seguintes valores: Lfv= 890 mm = 0,890 m Ltv = 1582 mm = 1,582 m Hgv = 649 mm = 0,649 m 33 3.2 DETERMINAÇÃO DOS COEFICIENTES DE RIGIDEZ Considerando as dificuldades de se obter os dados de rigidez para as molas de suspensão do veículo através do fabricante, optou-se pela determinação experimental dos mesmos. Foram selecionadas cinco pessoas pesando pelo menos 70 kg cada um, as quais foram pesadas na balança de uma farmácia, obtendo-se os seguintes resultados: A) 88,7 kg = 870 N B) 85,0 kg = 834 N C) 94,1 kg = 923 N D) 72,8 kg = 714 N E) 89,7 kg = 880 N Os pneus do veículo foram calibrados todos com 29 psig e levado a uma superfície plana e horizontal, e preparado como segue: o centro de rotação de cada roda foi determinado, levantando-a com um macaco hidráulico, girando-a livremente e traçando com uma caneta fixa que produziu um pequeno círculo. O centro do círculo corresponde ao centro da roda. Na carroçaria do veículo, foram marcados pontos acima de cada roda do veículo (Figuras 7 e 8). Os pontos foram identificados como seguem: • Roda dianteira esquerda 1 • Roda dianteira direita 2 • Roda traseira esquerda 3 • Roda traseira direita 4 • Ponto no chassi sobre roda 1: 5 • Ponto no chassi sobre roda 2: 6 • Ponto no chassi sobre roda 3: 7 • Ponto no chassi sobre roda 4: 8 34 Figura 7 - Numeração dos pontos de referência no veículo – Vista em Elevação Figura 8 - Numeração dos pontos de referência no veículo – Vista em Planta Com o veículo vazio, mediram-se as alturas dos centros de cada roda e das marcas dos chassis em relação ao solo, utilizando-se uma escala de aço. Em seguida, os voluntários foram dispostos no veículo conforme Figuras 9 e 10. Figura 9 - Carregamento do veículo com passageiros – Vista em Elevação 35 Figura 10 - Carregamento do veículo com passageiros – Vista em Planta Foram determinadas as posições médias dos Centros de Gravidade dos ocupantes, em relação ao eixo dianteiro, obtendo-se as dimensões Lfa,b e Lfc,d,e como na Figura. 11, bem como as alturas Ha,b e Hc,d,e dos CG dos mesmos em relação ao solo. Para a posição do CG dos ocupantes, considerou-se o centro do tronco na altura do umbigo. Com esses dados, foram calculadas as cargas de apoio das rodas nos eixos dianteiro e traseiro, conforme Figura 11. Em seguida, foram medidas as alturas dos pontos 1 a 8 em relação ao solo, para obter a deflexão de cada ponto com as cargas adicionais, e determinar os coeficientes de rigidez das molas de suspensão e dos pneus Figura 11 - Determinação das reações de apoio dianteiro e traseiro com o veículo carregado Através da somatória de momentos em relação aos pontos de apoio das rodas foram obtidas as reações de apoio Wkf e Wkt que foram consideradas igualmente distribuídas às rodas da direita e esquerda , frontais e traseiras, respectivamente. Além disso, calculou-se a posição do Centro de Gravidade do 36 veículo carregado, também pelo método de somatória de momentos (Taborek, 1957), obtendo-se os valores: Peso total do veículo carregado: Wk = 14028 N Distância do CG do veículo carregado ao eixo dianteiro: Lfk = 1,135 m 3.2.1. Determinação do Coeficiente de Rigidez dos Pneus Na Figura 12, é representada uma roda genérica, inicialmente submetida a uma carga Wfv1 relativa ao veículo vazio, e nessa situação a medida do centro até o solo é fvi, sendo i = 1 a 4, correspondendo à numeração atribuída às rodas no veículo. Em seguida, a mesma roda é submetida a uma carga Wfk1, relativa ao veículo carregado, sofrendo uma deflexão, e a medida da distância do centro da roda ao solo passa a ser fki. Portanto, a deflexão da roda com o veículo carregado deve ser calculada com: dfi = fv i – fk i (4) Figura 12 – Esquema para determinação do coeficiente de rigidez do pneu Os valores das deflexões nos pontos determinados sobre o veículo estão indicados na Tabela 1. 37 Tabela 1 – Deflexões das suspensões e rodas do Clio (m) Centros das rodas Marcas no veículo frontal traseira frontal traseira Pontos de referência 1 2 3 4 5 6 7 8 Veículo vazio fv1 fv2 fv3 fv4 fv5 fv6 fv7 fv8 0,265 0,265 0,273 0,270 0,622 0,618 0,620 0,626 Veículo carregado com 5 voluntários, conforme Figuras 9 e 10 fk1 fk2 fk3 fk4 fk5 fk6 fk7 fk8 0,257 0,261 0,260 0,262 0,590 0,587 0,535 0,545 Deflexão verificada (m) df1 df2 df3 df4 df5 df6 df7 df8 0,008 0,004 0,013 0,008 0,032 0,031 0,085 0,081 Nota-se na Tabela 1 que as deflexões nos lados esquerdos do veículo, tanto à frente como atrás, nas rodas e nas marcas da carroçaria do veículo, são maiores do que no lado direito. A diferença na distribuição de peso dos voluntários para o lado direito e lado esquerdo, em torno de 4,4%, não explicaria tal comportamento, principalmente para os pneumáticos. A fim de permitir a obtenção de um valor aproximado dos coeficientes de rigidez para as simulações, considerou-se a menor deflexão entre os lados direito e esquerdo, e distribuição de peso igualmente distribuído entre os lados também. Assim, para a deflexão dos pneus foi adotado o seguinte valor: • Rodas dianteiras: dff = df2 = 0,04 m • Rodas traseiras: dft = df4 = 0,08 m E para as cargas em cada roda, foram adotados os seguintes valores: Com o veículo vazio: • Carga dianteira média por roda: 1 WfvWfv 2 = (5) Wfv1 = 320 kgf = 3139 N • Carga traseira média por roda: 1 WtvWtv 2 = (6) Wtv1 = 180 kgf = 1766 N 38 Com o veículo carregado: • Carga dianteira média por roda: 1 WfkWfk 2 = (7) Wfk1 = 386,8 kgf = 3795 N • Carga traseira média por roda: 1 WtkWtk 2 = (8) Wtk1 = 328,4 kgf = 3221,6 N O coeficiente de rigidez de cada pneu é determinado como segue: Pneu dianteiro: 1 1Wfk WfvKPf dff − = (9) KPf = 164000 N/m Pneu traseiro: 1 1Wtk WtvKPt dft − = (10) KPt = 181950 N/m Os resultados apresentaram diferença de 6,76 %, tendo-se adotado o valor médio: KP = 173000 N/m 3.2.2 Determinação do Coeficiente de Rigidez das Molas de Suspensão Para a determinação do coeficiente de rigidez das molas de suspensão dianteira e traseira do veículo, deve-se inicialmente descontar a parcela de deformação do pneumático, permanecendo somente a deformação da suspensão. Na Figura 13 pode-se observar que quando o veículo é carregado, a carga inicial Wfv1 passa para Wfk1, a massa suspensa sofre um deslocamento df6, mas o pneu também sofre um achatamento df2, como foi mostrado na Figura 11. 39 Portanto, a diferença entre os dois valores corresponde à deformação somente da mola de suspensão. Figura 13 – Esquema para determinação do coeficiente de rigidez da mola de suspensão principal 3.2.2.1 Coeficiente de Rigidez das Molas de Suspensão Frontal Considerando-se o acima exposto, pode-se determinar a deflexão da suspensão principal dianteira do veículo, pela seguinte expressão: 6 2 6 2df df df− = − (11) ∆df6-2 = 0,027 m O aumento de carregamento é dado por: 1 1Wf Wfk Wfv= − (12) ∆Wf = 656 N Calcula-se o coeficiente de rigidez aparente da mola de suspensão dianteira, no plano da roda pela expressão: app 6 2 WfKF df − = (13) KFapp = 24296 N/m Como a mola principal não está montada diretamente sobre a roda, considera-se a relação entre os braços de suspensão e o ponto de fixação do amortecedor ao braço conforme mostrado na Figura 14. 40 Figura 14 - Montagem esquemática da suspensão da roda dianteira Assim, calcula-se o coeficiente de rigidez efetivo da mola principal com a expressão: app 0,3225KF KF 0,180 = ⋅ (14) KF = 43487 N/m Adota-se o valor: KF = 43500 N/m. 3.2.2.2 Coeficiente de Rigidez das Molas de Suspensão Traseira Para determinar o coeficiente de rigidez das molas de torção das rodas traseiras, considera-se a Figura 16, onde a roda (4) é representada inicialmente carregada com Wfv1, correspondendo ao veículo vazio. A distância do centro da roda ao solo é fv4, e a do centro de articulação do braço de suspensão ao solo é fv8. Em seguida, o mesmo conjunto está submetido à carga Wfk1, correspondendo ao veículo carregado. As distâncias correspondentes passam a ser fk4 e fk8. Tratando-se de mola de torção, o cálculo do coeficiente de rigidez do sistema de suspensão traseiro deve considerar o momento de torção aplicado e a deformação angular do braço. 41 Figura 15 - Esquema da suspensão traseira, veículo vazio e carregado O deslocamento vertical da suspensão do veículo, de vazio a carregado, é obtido considerando-se a diferença entre o deslocamento vertical do veículo e o da roda, similar ao adotado para a suspensão dianteira, como segue: 8 4 8 4df df df− = − (14) ∆df8-4 = 0,073 m O coeficiente de rigidez a torção foi calculada segundo a expressão abaixo (Roark, 1975), efetuando-se os necessários ajustes de unidades. TKT = θ (15) Wt L4KT ⋅ = θ (16) Considera-se que o valor de L4 é constante para pequenos valores de ∆θ: Temos que: 1 1Wt Wtk Wtv= − (17) ∆Wt = 1455,6 N 8 4dftg L4 −θ = (18) tg∆ θ = 0,1921 Então, ∆θ = 10,8741 grau = 0,1898 rad 42 Substituindo-se os valores obtidos bem como o valor de L4 = 0,38 m na equação (16), temos: KT = 2914 N.m/rad 3.2.2.3 Dados dos Amortecedores Na falta de dados do fabricante, foram realizadas algumas simulações e adotados valores adequados de coeficientes de amortecimento para a pior situação, seguido de fixação desses valores e estudo de comportamento nas demais situações. Assim, foram fixados os seguintes valores: Amortecedores dianteiros: CD = 1800 N.s/m Amortecedores traseiros: CT = 1100 N.s/m Os cursos dos amortecedores são de 0,180 m, para os dianteiros e 0,200 m para os traseiros, conforme medições feitas nos mesmos. 43 4. SIMULAÇÕES 4.1 PROCESSAMENTO NO PROGRAMA FORTRAN O programa AUTOLEV processa o arquivo Clio63.atl e emite um programa em FORTRAN, denominado Clio63.FOR e cria alguns arquivos auxiliares, entre eles o arquivo de entrada, denominado clio63.in, onde se digitam os dados na seqüência prevista. A listagem do programa em FORTRAN encontra-se no Apêndice 2. Ao executar o programa, define-se o tempo de tráfego, os intervalos de tempo a considerar e os passos de processamento, este último sendo cerca de 1/10 dos intervalos de tempo. O programa FORTRAN tem como resultados os deslocamentos previstos nos sete graus de liberdade, distribuídos em diversos arquivos em forma de tabelas. 4.2 PLOTAGEM DAS CURVAS Os resultados são plotados em curvas utilizando-se um programa adequado, que permite a seleção prévia e agrupamento dos dados por tipo de simulação, etc., permitindo-se a análise dos gráficos e introduzindo-se as alterações nos dados de entrada para novas simulações. 4.3 CONDIÇÕES DE TRÁFEGO CONSIDERADAS NAS SIMULAÇÕES Selecionaram-se as seguintes rugosidades, sendo uma com amplitude relativamente alta e a outra com rugosidade menor e passo menor: • Estradas precárias, clima severo – amplitude de 0,050 m com passos de 0,75 m; • Estradas precárias– amplitude de 0,025 m e passos de 0,30 m; Observa-se que outras ondulações encontradas em literaturas, consideradas longas com mais de 6 metros de passo não estão incluídas no presente estudo por estarem muito longe do objetivo inicial. 44 Consideram-se as seguintes velocidades de tráfego: 30 km/h, 60 km/h, 80 km/h e120 km/h (8,33 m/s, 16,67 m/s, 22,22 m/s e 33,33 m/s respectivamente). Estudou-se situações de carregamento do veículo com uma pessoa e com cinco pessoas, alterando-se nesse caso a posição do centro de gravidade S do corpo do veículo. Realizou-se as simulações utilizando-se uma tabela de DADOS DE ENTRADA, onde indicando-se todos os dados na seqüência de entrada. Para facilitar a visualização dos dados que serão alterados a cada simulação, apenas a primeira coluna é completa, correspondendo ao primeiro ensaio, e indica-se somente os valores que devem ser alterados a cada simulação em cada coluna subseqüente. As simulações estão separadas em duas séries: Simulações Série 5: para seqüências de 5 (cinco) passos de onda: Simulações Série 6: para seqüências de 20 (vinte) passos de onda: 4.4 RESULTADOS DE SIMULAÇÕES No item seguinte, são analisados os principais resultados, com o objetivo de verificar as situações de oscilações muito próximas ou acima dos valores limites de curso, esforços superiores ao peso do veículo em relação ao ponto analisado, tendências de ampliação dos picos de oscilação e de esforços ao longo do trajeto nas ondulações, e outros efeitos observados. As simulações não comentadas neste item estão apresentadas na forma gráfica no Apêndice 5. 45 4.5 ANÁLISE DOS RESULTADOS Para auxiliar as análises, calculam-se os intervalos de tempo gastos para a passagem das rodas pelas rugosidades: • Tempo inicial t1, correspondente ao intervalo de tempo entre o primeiro contato das rodas dianteiras e o primeiro contato das rodas traseiras com a rugosidade da pista; • Tempo t2, correspondente ao intervalo de tempo em que o veículo inteiro está inteiramente dentro do trecho ondulado; • Tempo t3: correspondente ao intervalo de tempo entre o primeiro contato das rodas dianteiras e o primeiro contato das rodas traseiras com o trecho plano da pista após o trecho rugoso; No caso de pista ondulada com 5 passos de 0,75 m, temos uma extensão de 3,75 m e como a distância entre eixos é de 2,472 m, também devemos considerar um comprimento inicial de 2,472 m na entrada e igual comprimento na saída. O tempo gasto para esses percursos está representado na Tabela 2 a seguir: Tabela 2 – Tempos de percurso para o trecho com 5 ondulações Velocidade m/s Tempo t1 (s) Tempo t2 (s) Tempo t3 (s) Tempo total 8,33 0,30 s 0,45 s 0,30 s 1,05 s 16,67 0,15 s 0,23 s 0,15 s 0,53 s 22,22 0,11 s 0,17 s 0,11 s 0,39 s 33,33 0,07 s 0,11 s 0,07 s 0,25 s No caso de pista ondulada com 20 passos de 0,75 m, temos uma extensão de 15 m e como a distância entre eixos é de 2,472 m, também devemos considerar um comprimento inicial de 2,472 m na entrada e igual comprimento na saída. O tempo gasto para esses percursos está representado na Tabela 3 a seguir: 46 Tabela 3 – Tempo de percurso para trecho com 20 ondulações Velocidade m/s Tempo t1 (s) Tempo t2 (s) Tempo t3 (s) Tempo total 8,33 0,30 s 1,8 s 0,30 s 2,4 s 16,67 0,15 s 0,9 s 0,15 s 1,2 s 22,22 0,11 s 0,68 s 0,11 s 0,9 s 33,33 0,07 s 0,45 s 0,07 s 0,59 s Verifica-se nas simulações com o veículo vazio e carregado que o tráfego pelas rugosidades a 8,33 m/s (30 km/h) gera oscilações muito fortes, sendo que nas velocidades mais altas as oscilações são significativamente menores. Na Figura 16 - Gráfico DYA10-501, verifica-se que as oscilações verticais da roda dianteira à velocidade de 8,33 m/s em ondulações com amplitude de 0,05 m e passos de 0,75 m e seqüência de 5 passos, variam de 0,15 a -0,15 m, totalizando 0,3 m de extensão. Como o amortecedor está fixado em um ponto intermediário do braço, o curso correspondente é de 0,167 m, valor muito próximo ao curso total do amortecedor que é de 0,180 m, podendo atingir um dos limites de curso. 0 1 2 3 -0,15 -0,10 -0,05 0,00 0,05 0,10 0,15 0,20 GRAFICO DYA10-501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupanteos ci la çã o ve rti ca l d a ro da d ia nt ei ra (m ) tempo (s) DYA10-501 v= 8,33 m/s DYA10-502 v= 16,67 m/s DYA10-503 v= 22,22 m/s DYA10-504 v= 33,33 m/s Figura 16 - Gráfico DYA10- 501 A 504 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 47 Estando o veículo lotado com cinco ocupantes conforme Figura 17 - Gráficos DYA10-505 A 508, observa-se maiores oscilações do braço em todas as velocidades, principalmente a 8,33 m/s e 16,67 m/s, amortecendo após a saída do trecho ondulado a 1,05 s. 0 1 2 -0,02 0,00 0,02 0,04 0,06 GRAFICO DYA10-505 A 508 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo cinco ocupantes O sc ila çã o ve rti ca l d a ro da d ia nt ei ra (m ) Tempo (s) DYA10-505 v= 8,33 m/s DYA10-506 v= 16,67 m/s DYA10-507 v= 22,22 m/s DYA10-508 v= 33,33 m/s Figura 17 - Gráfico DYA10- 505 A 508 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos, veículo com cinco pessoas. Por outro lado, pode-se verificar na Figura 18 - Gráfico DYA10-601 A 604, que trafegando a 8,33 m/s em trecho com as mesmas características e seqüência de 20 passos de ondulação, os picos chegam a 0,22 e -0,14 m, com extensão de 0,36 m, sendo 0,2 m para o amortecedor, o que significa que os limites de curso da suspensão foram atingidos. Na Figura 19 - Gráfico DYA11-501 a 504, as oscilações da roda traseira para a velocidade de 8,33 m/s nas mesmas ondulações variam entre 0,25 e -0,2 m, dando um total de 0,45 m, correspondendo a 0,25 m no amortecedor, significando que há impacto no final de curso. 48 0 1 2 3 4 -0,15 -0,10 -0,05 0,00 0,05 0,10 0,15 0,20 0,25 GRAFICO DYA10 - 601 A 604 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante O sc ila çã o ve rti ca l d as ro da s di an te ira s (m ) Tempo (s) DYA10-601 v= 8,33 m/s DYA10-602 v= 16,67 m/s DYA10-603 v= 22,22 m/s DYA10-604 v= 33,33 m/s Figura 18 - Gráfico DYA10- 601 A 604 – Oscilação vertical das rodas dianteiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 0 1 2 3 -0,2 -0,1 0,0 0,1 0,2 0,3 GRAFICO DYA11 - 501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante O sc ila çã o ve rti ca l r od a tra se ira e sq ue rd a (m ) Tempo (s) DYA11-501 v= 8,33 m/s DYA11-502 v= 16,67 m/s DYA11-503 v= 22,22 m/s DYA11-504 v= 33,33 m/s Figura 19 - Gráfico DYA11- 501 A 504 – Oscilação vertical das rodas traseiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 49 Na Figura 20 - Gráfico DYA11-601 A 604, com 20 passos de ondulação, há tendência de amplificação das oscilações, chegando a 0,42 e -0,33m, ou seja, total de 0,75 m, resultando em 0,33 m para o amortecedor, sendo superior a 0,20 m, que é totalmente inaceitável. Pode-se notar também que as amplitudes de oscilação da roda traseira tende a aumentar, desde 0,3 s quando inicia a passagem pelo trecho ondulado, e o amortecimento somente ocorre a partir de 2,1 s quando a roda traseira entra no trecho plano. A freqüência de oscilação está entre 6,5 e 7 Hz. 0 1 2 3 4 -0,3 -0,2 -0,1 0,0 0,1 0,2 0,3 0,4 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante GRAFICO DYA11 - 601 A 604 O sc ila çã o ve rti ca l d as ro da s tra se ira s (m ) Tempo (s) DYA11-601 v= 8,33 m/s DYA11-602 v= 16,67 m/s DYA11-603 v= 22,22 m/s DYA11-604 v= 33,33 m/s Figura 20 - Gráfico DYA11- 601 A 604 – Oscilação vertical das rodas traseiras em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. Estando o veículo apenas com o motorista e trafegando na ondulação maior, o primeiro contato das rodas dianteiras com a rugosidade gera uma força de compressão no ponto de ligação da suspensão com o monobloco, de quase 15000 N em uma roda (Figura 21 - Gráfico FA1-501a504), significando que o veículo perde o contato com o solo devido ao impacto, uma vez que na suspensão dianteira a carga por roda com uma pessoa não ultrapassa 3400 N 50 0 1 2 3 4 -15000 -10000 -5000 0 5000 GRAFICO FA1 - 501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante Fo rç a na c on ex ão s up er i0 or d o am or te ce do r d ia nt ei ro (N ) Tempo (s) FA1-501 v= 8,33 m/s FA1-502 v= 16,67 m/s FA1-503 v= 22,22 m/s FA1-504 v= 33,33 m/s Figura 21 - Gráfico FA1- 501 A 504 – Força na conexão superior do amortecedor dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. À mesma velocidade e em trecho com 20 passos de ondulação, pode-se notar na Figura 22 - Gráfico FA1-601 A 604, que há tendência de surgirem picos de carga a intervalos de 0,5 s (2Hz) com tendência a crescimento, chegando a atingir 8000 N e -14000 N até que a roda traseira entra no trecho plano a 2,4 s e ainda assim a oscilação continua até que, depois de 0,5 s da última rugosidade a oscilação começa a diminuir. Para as demais velocidades consideradas, entretanto, todos os picos observados são menores, apresentando tendência de amortecimento mesmo no trecho em que todo o veículo está percorrendo as ondulações. 51 0 1 2 3 4 5 -15000 -10000 -5000 0 5000 10000 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante GRAFICO FA1 - 601 A 604 Fo rç a na c on ex ão s up er io r d o am or te ce do r d ia nt ei ro (N ) Tempo (s) FA1-601 v= 8,33 m/s FA1-602 v= 16,67 m/s FA1-603 v= 22,22 m/s FA1-604 v= 33,33 m/s Figura 22 - Gráfico FA1- 601 A 604 – Força na conexão superior do amortecedor traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos Situação semelhante ocorre no ponto de conexão do amortecedor da roda traseira. Percorrendo trecho com 5 ondulações, surgem picos de 10000 N em compressão e 8000 N de tração enquanto o veículo está totalmente dentro da região rugosa (Figura 23 - Gráfico FA3-501 A 504). Ao sair a 1,05 s o amortecimento ocorre em praticamente 0,5 s. Outrossim, na pista com 20 ondulações, na Figura 24 - Gráfico FA3-601 A 604, nota-se uma série de picos de compressão com a freqüência em torno de 9 a 15 Hz chegando a atingir 10000 N, e dois picos de tração atingindo 12000 N e 17000 N, com tendência de alta. 52 0 1 2 -10000 -8000 -6000 -4000 -2000 0 2000 4000 6000 8000 GRAFICO FA3 - 501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante Fo rç a na c on ex ão s up er io r d o am or te ce do r t ra se iro (N ) Tempo (s) FA3-501 v= 8,33 m/s FA3-502 v= 16,67 m/s FA3-503 v= 22,22 m/s FA3-504 v= 33,33 m/s Figura 23 - Gráfico FA3- 501 A 504 – Força na conexão superior do amortecedor dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 0 1 2 3 4 -10000 -5000 0 5000 10000 15000 GRAFICO FA3-601 A 604 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante Fo rç a na c on ex ão s up er io r d o am or te ce do r t ra se iro (N ) tempo (s) FA3-601 v= 8,33 m/s FA3-602 v= 16,67 m/s FA3-603 v= 22,22 m/s FA3-604 v= 33,33 m/s Figura 24 - Gráfico FA3- 601 A 604 – Força na conexão superior do amortecedor traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 53 Na suspensão traseira, o torque da mola de torção sobre o braço apresenta picos de 6600 N.m a -2500 N.m com o veículo totalmente dentro do trecho rugoso com 5 ondulações (Figura 25 - Gráfico TORQ-501 A 504). Considerando o braço da suspensão de 0,38 m, a carga correspondente na roda é de 17370 N e 6580 N respectivamente, muito superiores à carga em torno de 2000 N de cada roda traseira com o veículo transportando apenas uma pessoa. 0 1 2 3 -2000 0 2000 4000 6000 GRAFICO TORQ-501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante To rq ue n o br aç o de s us pe ns ão tr as ei ro (N .m ) tempo (s) TORQ-501 v= 8,33 m/s TORQ-502 v= 16,67 m/s TORQ-503 v= 22,22 m/s TORQ-504 v= 33,33 m/s Figura 25 - Gráfico TORQ-501 A 504 – Torque no braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. Quando percorre trecho com 20 ondulações (Figura 26 - Gráfico TORQ- 601 A 604), nota-se a tendência de aumento dos valores a cada pico, e o amortecimento ocorre imediatamente após a saída do trecho rugoso, a 2,4 segundos. 54 0 1 2 3 4 -4000 -2000 0 2000 4000 6000 8000 GRAFICO TORQ-601 A 604 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante To rq ue n o br aç o de s us pe ns ão tr as ei ro (N .m ) Tempo (s) TORQ-601 v= 8,33 m/s TORQ-602 v= 16,67 m/s TORQ-603 v= 22,22 m/s TORQ-604 v= 33,33 m/s Figura 26 - Gráfico TORQ- 601 A 604 – Torque no braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. A oscilação do braço de suspensão dianteiro apresentou uma particularidade (Figura 27 – Gráfico Q4-501 A 504). Após a passagem pelo trecho de amplitude 0,05 m e passo de 0,75 m a 16,67 m/s, ela se estabilizou em torno de 0,75 graus, diferente das outras curvas onde a estabilização ocorreu a cerca de 2,3 graus. Trafegando a 8,33 m/s, a rotação apresenta uma ampla oscilação logo após a saída da roda traseira das ondulações, antes da estabilização. A rotação do braço em 3 graus significa cerca de 10 mm de curso no amortecedor, que está a 0,180 m do ponto de articulação, podendo ser considerado desprezível. Como se trata de uma simulação virtual, esse comportamento pode ser atribuído a uma não-linearidade nas equações de movimento. Comportamento similar apresenta o braço de suspensão dianteiro na passagem pela seqüência de 20 passos (Figura 28 - Gráfico Q4-601 a 604), a velocidade de 8,33 m/s, após a qual estabiliza-se a -4,5 graus, enquanto que as demais curvas estabilizam-se a 2,3 graus. Esse comportamento se repete em diversas simulações, mesmo em computadores diferentes utilizados. 55 0 1 2 3 4 -0,5 0,0 0,5 1,0 1,5 2,0 2,5 3,0 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante GRAFICO Q4-501 A 504 R ot aç ão d o br aç o di an te iro (g ra u) Tempo (s) Q4-501 v= 8,33 m/s Q4-502 v= 16,67 m/s Q4-503 v= 22,22 m/s Q4-504 v= 33,33 m/s Figura 27 - Gráfico Q4- 501 A 504 – Rotação do braço de suspensão dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 0 1 2 3 4 -6 -4 -2 0 2 4 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante GRAFICO Q4-601 A 604 R ot aç ão d o br aç o de s us pe ns ão d ia nt ei ra (g ra u) Tempo (s) Q4-601 v= 8,33 m/s Q4-602 v= 16,6,7 m/s Q4-603 v= 22,22 m/s Q4-604 v= 33,33 m/s Figura 28 - Gráfico Q4-601 A 604 – Rotação do braço de suspensão dianteiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. Trafegando a 8,33 m/s pelas seqüências de 5 e de 20 passos, (Figura 29- Gráfico Q6-501 A 504 e Figura 30 – Gráfico Q6-601 A 604), observa-se também no braço de suspensão traseira, a tendência de aumento das oscilações, que se atenuam quando as rodas traseiras deixam o trecho ondulado, no primeiro caso a 56 1,05 s e no segundo a 2,4 s. Por outro lado, no braço de suspensão traseira, a rotação de 2 graus representa para o amortecedor localizado a 0,165 m do ponto de articulação um curso de 6 mm, que pode ser considerado desprezível. 0 1 2 3 -2,5 -2,0 -1,5 -1,0 -0,5 0,0 0,5 1,0 GRAFICO Q6-501 A 504 amplitude = 0,050 m passo = 0,750 m no. passos = 5 Veiculo um ocupante R ot aç ão d o br aç o de s us pe ns ão tr as ei ro (g ra u) Tempo (s) Q6-501 v= 8,33 m/s Q6-502 v= 16,67 m/s Q6-503 v= 22,22 m/s Q6-504 v= 33,33 m/s Figura 29 - Gráfico Q6-501 A 504 – Rotação do braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 5 passos. 0 1 2 3 4 -3 -2 -1 0 1 2 amplitude = 0,050 m passo = 0,750 m no. passos = 20 Veiculo um ocupante GRAFICO Q4-601 A 604 R ot aç ão d o br aç o de s us pe ns ão tr as ei ra (g ra u) tempo (s) Q6-601 v= 8,33 m/s Q6-602 v= 16,67 m/s Q6-603 v= 22,22 m/s Q6-604 v= 33,33 m/s Figura 30 - Gráfico Q6-601 A 604 – Rotação do braço de suspensão traseiro em trecho com ondulações de amplitude 0,05 m e passo de 0,75 m, seqüência de 20 passos. 57 Analisando todos os gráficos, verifica-se que na passagem do veículo pelas ondulações a baixa velocidade, em torno de 8,33 m/s o veículo sofre oscilações e esforços de grande amplitude. Esse comportamento ocorre com o veículo quase vazio, estando com apenas um ocupante, e também estando cheio com cinco pessoas. Como foi verificado que há ocorrência de oscilações superiores ao curso do amortecedor ou mesmo de valores superiores à carga correspondente em cada roda, a possibilidade de o veículo perder o contato com a pista é concreta. Por outro lado, a velocidades superiores a 16,67 m/s o veículo apresenta comportamento bem melhor, com ausência de grandes oscilações, e mesmo em trechos longos com 20 passos, não mostra tendência de oscilações crescentes. Conclui-se que a freqüência de excitação a 8,33 m/s aproxima-se da freqüência natural da suspensão. As oscilações de grande amplitude e repetitivas levam à redução das características dos amortecedores devido à solicitação intensa, que aumenta a temperatura do fluido e reduz a viscosidade, além da dissolução do gás ao fluido em alguns casos. 58 5 CONCLUSÕES A análise dos resultados obtidos nos ensaios leva a concluir que mesmo um veículo com suspensão passiva e projeto de suspensão bem equilibrada, é de se esperar que, em certas condições de tráfego por pavimentos que apresentem rugosidade muito regular pode apresentar comportamento de ressonância. De um modo geral, a passagem por trechos ondulados a velocidades relativamente baixas tende a excitar o veículo a freqüências muito próximas da freqüência natural da suspensão, tornando a condução bastante desconfortável, além da possibilidade de haver a perda de contato das rodas com o solo. No veículo de pequeno porte como é o caso em estudo, verifica-se que a freqüência natural da massa suspensa do veículo está em torno de 0,8 a 1 Hz, e a massa não suspensa apresenta freqüência natural da ordem de 3,5 a 4 Hz. Podem ser apresentadas as seguintes sugestões para melhorar as condições de tráfego do veículo em estudo: • Redimensionar os conjuntos de mola e amortecedor dianteiro e traseiro, para aumentar a rigidez e em paralelo o coeficiente de amortecimento de choques. O conforto para os ocupantes do veículo diminui, enquanto aumenta a segurança do tráfego. Essa solução também deve levar em conta o eventual reforço do monobloco para suportar cargas maiores que a mola irá submeter ao seu ponto de apoio, conforme Sharp (Sharp, 1986); • Modificar a relação entre o comprimento do braço e a posição de montagem da mola e do amortecedor no braço, para aumentar a eficiência. Evidentemente, isso significa também modificar o comprimento da mola e o curso do amortecedor, com a conseqüente necessidade de redimensionar o espaço necessário para a instalação; • A utilização de amortecedores com ajuste de amortecimento por meio eletrônico ou pneumático, embora não seja a proposta inicial deste 59 trabalho, poderia ser estudado caso for verificado que o custo de implantação estaria acessível para as faixas de veículos compactos. Embora os gráficos indiquem que o tráfego por pistas com ondulações a velocidades altas produz em geral menores oscilações, deve-se considerar que nessas condições, a dirigibilidade do veículo torna-se seriamente comprometida, com dificuldades em realizar curvas de raios relativamente pequenos, pois os contatos das rodas com a pista são reduzidos. Com o objetivo de contribuir para outras pesquisas no campo de dinâmica de veículos, apresenta-se a seguir algumas sugestões de assuntos, incluindo variações de abordagem, vem como estudos de soluções. No presente trabalho considerou-se que o veículo trafega em linha reta pelas ondulações da pista, ou seja, as rodas da direita e da esquerda são submetidas às mesmas condições de rugosidade. Com isso, o momento de rolamento q2 permaneceu nulo. Sendo assim, seria interessante introduzir defasagem na passagem das rodas da direita e da esquerda, simulando rugosidades diagonais e estudar o comportamento do veículo nessas condições. Outra abordagem interessante seria a definição de outros graus de liberdade, como o assento do motorista e de pelo menos de um dos passageiros do banco traseiro, procurando-se obter as oscilações sobre o ocupante nos sentidos vertical, rotação segundo o eixo transversal, ou seja a arfagem, bem como os esforços verticais atuantes. Embora não tenha sido o objetivo do presente trabalho, o estudo de soluções econômicas para suspensões ativas ou semi-ativas seria muito interessante. Nesse ponto de vista, amortecedores com ajuste de amortecimento por meio eletrônico ou pneumático com acionamento semi-automático ou mesmo manual poderia ser desenvolvido. O desenvolvimento de molas com rigidez variável controlável é desejável, devendo ser objeto de desenvolvimento. Há uma solução com mola pneumática, cuja rigidez é aumentada por injeção de ar comprimido, ou uma combinação de mola pneumática com mola de aço onde a mola pneumática aumenta a sua 60 parcela de rigidez por meio de injeção de ar comprimido. (Sharp, 1986). O custo desta solução poderia ser reduzido com controle manual do processo. No presente trabalho foram selecionados dois tipos de rugosidade, sendo que o segundo, com amplitude menor e passo também menor não apresentou oscilações importantes, mas estudos com outros tipos de rugosidades devem acrescentar informações interessantes para o estudo de suspensões passivas. 61 REFERÊNCIAS BASTOW, D., HOWARD, G. - Car Suspension and Handling, Third Edition, Pentech Press Ltd., London and Society of Automobile Engineers, Inc., 1993; BIANCHI,G., SCHIEHLEN,W.O., Editors, Dynamics of Multibody Systems, Italy, Springer-Verlag, 1985; DOBROVOLSKI et al, Elementos de Máquinas – pg. 174 – Editorial Mir – Moscú, 1980 GILLESPIE, T.D., Fundamentals os Vehicle Dynamics, Society of Automobile Engineers, Inc. 1992; KANE, T.R., Dynamics of Nonholonomic Systems, Journal of Applied Mechanics, vol. 28, n.4, pag. 574-578, 1961; KANE, T.R., LEVINSON, D.A., Multibody Dynamics, Journal of Applied Mechanics, vol. 50, pag. 1071-1078, 1978; KANE, T.R., LEVINSON, D.A., Formulation of Equation for Complex Spacecraft, Journal of Guidance and Control, vol. 3, n.2, pag. 99-112, 1980; KANE, T.R., LEVINSON, D.A., The Use of Kane’s Dynamical Equations in Robotics, The International Journal of Robotics Research, vol. 2, n. 3, pag. 3-21, 1983; KANE, T.R., LEVINSON, D.A., Dynamics : Theory and Applications, McGraw-Hill, USA, 1985; KÖRTUM, W., SCHIELEN,W., General Purpose Vehicle System Dynamics Software based on Multibody Formalisms Vehicle System Dynamics, vol. 14, pag. 229-263, 1985; KREUZER, E.J., SCHIEHLEN, W.O., Equations of Motion and Equations of Stres for Robots and Manipulators: An Application of the NEWEUL formalism, Proceedings of RoManSy’84 : The Fifth CISM-IFToMM Symposium : Theory and Practice of Robots and Manipulators, pag. 79-85, 1984; 62 MATHIAS, M.H. – Notas de Aula – Análise dos Movimentos da Suspensão Primária / Marcha Suave– Faculdade de Engenharia de Guaratinguetá – UNESP, 2004 NEWTON, K., STEEDS, W., GARRET, T.K., The Motor Vehicle, 11a. Edition, Butterworth-Heinemann, 1989; NIEMANN, G. – Elementos de Máquinas Volume I – Ed. Edgard Blücher, 1971; ROBSON, J. D. – Road Surface Description and Vehicle Response – Int. Jdesign, vol. I, no. 1, Printed in U.K., 1979; ROARK, R.J., YOUNG, W.C. - Formulas for Stress and Strain, Fifth edition – chap.9, Table 20, McGraw-Hill–Kogakusha, 1975 SCHIEHLEN, W.O., KREUZER, E.J., Simbolic Computerized Derivation of Equation of Motion, Dynamics of Multibody Systems Dynamics, Symposium Munich/Germany, Editor K. Magnus, pag. 290-305, 1977; SCHIEHLEN, W.O., KREUZER, E.J., Strength Estimations in Multibody Systems, Dynamics of Multibody Systems Dynamics, Symposium Udine/Italy, Editors G. Bianchi and W. Schiehlen, pag. 249-259, 1985; SHABANA, A.A., Automated Analysis of Constrained Systems of Rigid and Flexible Bodies, J. Vibration, Acoustics, Stress and Reability in Design, vol. 107, pag. 431-439, 1985; SHABANA, A.A., Dynamics of Inertia-Variant flexible Systems Using Experimentally Identified Parameters, J. Mechanisms, Transmissions, and Automation in Design, vol. 108, pag. 358-366, 1986. SHARP R.S. and HASSAN, S.A. – An Evaluation of Passive Automotive Suspension Systems with Variable Stiffness and Damping Parameters - Vehicle System Dynamics, 15 pp. 335-350, 1986 SHARP R.S. and CROLLA, D.A. – Road Vehicle Suspension System Design – A Review – Vehicle System Dynamics, 16 pp. 167-192, 1987 SHABANA, A.A., Dynamics of Multibody Systems, John Wiley & Sons, New York, USA, 1989 TABOREK, J.J. – Mechanics of Vehicles – Pentech Press – Design, 1957 63 TOMAZINI, J.E. – O Modelo Multicorpo Aplicado a um Manipulador Modelo Rígido e Flexível – Tese de Doutorado - Escola de Engenharia de São Carlos, 1996; WILLIAMS, R.A. – Automotive Active Suspensions – Part 1: Basic Principles – D00297 IMechE, 1997 WILLIAMS, R.A. – Automotive Active Suspensions – Part 2: Practical Considerations – D00397 IMechE, 1997 64 APÊNDICES APÊNDICE 1 – Listagem do programa AUTOLEV para a versão Clio6.3 APÊNDICE 2 – Listagem do programa FORTRAN para a versão Clio 6.3 APÊNDICE 3 – CLIO – DADOS DE ENTRADA – SÉRIE 5 – PARA QUANTIDADE DE ONDAS = 5 APÊNDICE 4 – CLIO – DADOS DE ENTRADA – SÉRIE 6 – PARA QUANTIDADE DE ONDAS = 20 APÊNDICE 5 – Gráficos dos ensaios. APÊNDICE1 – Listagem do programa AUTOLEV para a versão Clio6.3 !PROGRAMA: CLIO63 ! DOF(7) FRAMES(S,X,D,H,P,R) POINTS(O,A1,A2,A3,A4,A5,A6,A7,A8,A9,A10,A11,A12,A13,A14,A15,A16) MASS(S,MS,A9,M9,A10,M9,A11,M11,A12,M11,D,MD,H,MD,P,MP,R,MP) MASSLESS(O,A1,A2,A3,A4,A5,A6,A7,A8,A13,A14,A15,A16,X) INERTIA(S,IS1,IS2,IS3,0,0,0) INERTIA(D,ID,0,ID,0,0,0) INERTIA(H,ID,0,ID,0,0,0) INERTIA(P,0,IP,IP,0,0,0) INERTIA(R,0,IP,IP,0,0,0) CONST(DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L1,L2,L3,L4) CONST(KD,KT,KP,CD,CT) CONST(V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FATOR3,FATOR4) ! !EQUAÇÕES DIFERENCIAIS CINEMATICAS Q1'=U1 Q2'=U2 Q3'=U3 Q4'=U4 Q5'=U5 Q6'=U6 Q7'=U7 ! !MATRIZES DE TRANSFORMAÇÃO ! SIMPROT(N,X,2,Q2) SIMPROT(X,S,1,Q3) DIRCOS(N,S) SIMPROT(S,D,1,Q4) SIMPROT(S,H,1,Q5) SIMPROT(S,P,2,Q6) SIMPROT(S,R,2,Q7) DIRCOS(N,D) DIRCOS(S,D) DIRCOS(X,D) DIRCOS(N,H) DIRCOS(S,H) DIRCOS(X,H) DIRCOS(N,P) DIRCOS(S,P) DIRCOS(X,P) DIRCOS(N,R) DIRCOS(S,R) DIRCOS(X,R) ! !VETORES DE POSIÇÃO ! POSSTAR=Q1*N3 PSSTARA1=DX1*S1-DY1*S2-DZ1*S3 PSSTARA2=DX1*S1+DY1*S2-DZ1*S3 PSSTARA3=-DX3*S1-DY3*S2+DZ3*S3 PSSTARA4=-DX3*S1+DY3*S2+DZ3*S3 PSSTARA5=DX5*S1-DY5*S2+DZ5*S3 PSSTARA6=DX5*S1+DY5*S2+DZ5*S3 PSSTARA7=-DX7*S1-DY7*S2+DZ7*S3 PSSTARA8=-DX7*S1+DY7*S2+DZ7*S3 PA5A13=-L2*D2 PA5A9=-L1*D2 PA5DSTAR=-(L1/2)*D2 PA6A14=L2*H2 PA6A10=L1*H2 PA6HSTAR=(L1/2)*H2 PA7A15=-L3*P1 PA7A11=-L4*P1 PA7PSTAR=-(L4/2)*P1 PA8A16=-L3*R1 PA8A12=-L4*R1 PA8RSTAR=-(L4/2)*R1 ! !VETORES POSIÇÃO ANTES DE QUALQUER DEFORMAÇÃO ! VETA1=DX1*N1-DY1*N2-DZ1*N3 VETA2=DX1*N1+DY1*N2-DZ1*N3 VETA3=-DX3*N1-DY3*N2+DZ3*N3 VETA4=-DX3*N1+DY3*N2+DZ3*N3 VETA9=DX5*N1-(DY5+L1)*N2+DZ5*N3 VETA10=DX5*N1+(DY5+L1)*N2+DZ5*N3 VETA11=-(DX7+L4)*N1-DY7*N2+DZ7*N3 VETA12=-(DX7+L4)*N1+DY7*N2+DZ7*N3 VETA13=DX5*N1-(DY5+L2)*N2+DZ5*N3 VETA14=DX5*N1+(DY5+L2)*N2+DZ5*N3 VETA15=-(DX7+L3)*N1-DY7*N2+DZ7*N3 VETA16=-(DX7+L3)*N1+DY7*N2+DZ7*N3 ! !COMPONENTES DOS VETORES ACIMA NA DIREÇÃO N3 ! YYA1=DOT(VETA1,N3) YYA2=DOT(VETA2,N3) YYA3=DOT(VETA3,N3) YYA4=DOT(VETA4,N3) YYA9=DOT(VETA9,N3) YYA10=DOT(VETA10,N3) YYA11=DOT(VETA11,N3) YYA12=DOT(VETA12,N3) YYA13=DOT(VETA13,N3) YYA14=DOT(VETA14,N3) YYA15=DOT(VETA15,N3) YYA16=DOT(VETA16,N3) ! !VETORES DE POSIÇÃO DOS PONTOS EM RELAÇÃO AO SISTEMA NEWTONIANO ! POA1=ADD(POSSTAR,PSSTARA1) POA2=ADD(POSSTAR,PSSTARA2) POA3=ADD(POSSTAR,PSSTARA3) POA4=ADD(POSSTAR,PSSTARA4) POA5=ADD(POSSTAR,PSSTARA5) POA6=ADD(POSSTAR,PSSTARA6) POA7=ADD(POSSTAR,PSSTARA7) POA8=ADD(POSSTAR,PSSTARA8) ! !VELOCIDADES ! VSSTARN=V*N1+Q1'*N3 WSX=Q3'*S1 WXN=Q2'*N2 WSN=ADD(WSX,WXN) WDS=Q4'*D1 WDN=ADD(WDS,WSN) WHS=Q5'*H1 WHN=ADD(WHS,WSN) WPS=Q6'*P2 WPN=ADD(WPS,WSN) WRS=Q7'*R2 WRN=ADD(WRS,WSN) V2PTS(N,S,SSTAR,A1) V2PTS(N,S,SSTAR,A2) V2PTS(N,S,SSTAR,A3) V2PTS(N,S,SSTAR,A4) V2PTS(N,S,SSTAR,A5) V2PTS(N,S,SSTAR,A6) V2PTS(N,S,SSTAR,A7) V2PTS(N,S,SSTAR,A8) V2PTS(N,D,A5,A13) V2PTS(N,D,A5,A9) V2PTS(N,D,A5,DSTAR) V2PTS(N,H,A6,A14) V2PTS(N,H,A6,A10) V2PTS(N,H,A6,HSTAR) V2PTS(N,P,A7,A15) V2PTS(N,P,A7,A11) V2PTS(N,P,A7,PSTAR) V2PTS(N,R,A8,A16) V2PTS(N,R,A8,A12) V2PTS(N,R,A8,RSTAR) ! !ACELERACOES ! ASSTARN=DERIV(VSSTARN,T,N) ADSTARN=DERIV(VDSTARN,T,N) AHSTARN=DERIV(VHSTARN,T,N) APSTARN=DERIV(VPSTARN,T,N) ARSTARN=DERIV(VRSTARN,T,N) ALFSN=DERIV(WSN,T,N) ALFDN=DERIV(WDN,T,N) ALFHN=DERIV(WHN,T,N) ALFPN=DERIV(WPN,T,N) ALFRN=DERIV(WRN,T,N) AA9N=DERIV(VA9N,T,N) AA10N=DERIV(VA10N,T,N) AA11N=DERIV(VA11N,T,N) AA12N=DERIV(VA12N,T,N) ! !DESLOCAMENTO NOS PONTOS SUBMETIDOS A ESFORCOS ! YA1=DOT(POA1,N3) YA2=DOT(POA2,N3) YA3=DOT(POA3,N3) YA4=DOT(POA4,N3) POA13=ADD(POA5,PA5A13) YA13=DOT(POA13,N3) POA9=ADD(POA5,PA5A9) YA9=DOT(POA9,N3) POA14=ADD(POA6,PA6A14) YA14=DOT(POA14,N3) POA10=ADD(POA6,PA6A10) YA10=DOT(POA10,N3) POA15=ADD(POA7,PA7A15) YA15=DOT(POA15,N3) POA11=ADD(POA7,PA7A11) YA11=DOT(POA11,N3) POA16=ADD(POA8,PA8A16) YA16=DOT(POA16,N3) POA12=ADD(POA8,PA8A12) YA12=DOT(POA12,N3) ! !VELOCIDADES NOS PONTOS SUBMETIDOS A ESFORÇOS ! YA1P=DERIV(YA1,T) YA2P=DERIV(YA2,T) YA3P=DERIV(YA3,T) YA4P=DERIV(YA4,T) YA9P=DERIV(YA9,T) YA10P=DERIV(YA10,T) YA11P=DERIV(YA11,T) YA12P=DERIV(YA12,T) YA13P=DERIV(YA13,T) YA14P=DERIV(YA14,T) YA15P=DERIV(YA15,T) YA16P=DERIV(YA16,T) ! !DESLOCAMENTOS RELATIVOS ! DYA1=RIGHT(YA1)-RIGHT(YYA1) DYA2=RIGHT(YA2)-RIGHT(YYA2) DYA3=RIGHT(YA3)-RIGHT(YYA3) DYA4=RIGHT(YA4)-RIGHT(YYA4) DYA9=RIGHT(YA9)-RIGHT(YYA9) DYA10=RIGHT(YA10)-RIGHT(YYA10) DYA11=RIGHT(YA11)-RIGHT(YYA11) DYA12=RIGHT(YA12)-RIGHT(YYA12) DYA13=RIGHT(YA13)-RIGHT(YYA13) DYA14=RIGHT(YA14)-RIGHT(YYA14) DYA15=RIGHT(YA15)-RIGHT(YYA15) DYA16=RIGHT(YA16)-RIGHT(YYA16) ! ! ! FORCAS E TORQUES ! FORCE(SSTAR)=MS*G*N3 FORCE(A13/A1)=-KD*(RIGHT(DYA1)-RIGHT(DYA13))*N3-CD*(RIGHT(YA1P)-RIGHT(YA13P))*N3 FORCE(A14/A2)=-KD*(RIGHT(DYA2)-RIGHT(DYA14))*N3-CD*(RIGHT(YA2P)-RIGHT(YA14P))*N3 FORCE(A9)=-KP*(RIGHT(DYA9)-ZZ1)*N3+M9*G*N3 FORCE(A10)=-KP*(RIGHT(DYA10)-ZZ2)*N3+M9*G*N3 FORCE(A15/A3)=-CT*(RIGHT(YA3P)-RIGHT(YA15P))*N3 FORCE(A16/A4)=-CT*(RIGHT(YA4P)-RIGHT(YA16P))*N3 FORCE(A11)=-KP*(RIGHT(DYA11)-ZZ3)*N3+M11*G*N3 FORCE(A12)=-KP*(RIGHT(DYA12)-ZZ4)*N3+M11*G*N3 FORCE(DSTAR)=MD*G*N3 FORCE(HSTAR)=MD*G*N3 FORCE(PSTAR)=MP*G*N3 FORCE(RSTAR)=MP*G*N3 TORQUE(S/P)=-KT*Q6*P2 TORQUE(S/R)=-KT*Q7*R2 FA1=-KD*(RIGHT(DYA1)-RIGHT(DYA13))-CD*(RIGHT(YA1P)-RIGHT(YA13P)) FA2=-KD*(RIGHT(DYA2)-RIGHT(DYA14))-CD*(RIGHT(YA2P)-RIGHT(YA14P)) FA3=-CT*(RIGHT(YA3P)-RIGHT(YA15P)) FA4=-CT*(RIGHT(YA4P)-RIGHT(YA16P)) TORQ=DOT(-KT*Q6*P2,N2) ! ! FORCA ATIVA GENERALIZADA FR ! ! FORCA DE INERCIA GENERALIZADA FRSTAR ! ! EQUACOES DO MOVIMENTO KANE CONTROLS(ZZ1,ZZ2,ZZ3,ZZ4,DYA9,DYA10,DYA11,DYA12,FA1,FA2,FA3,FA4,TORQ) ! ! UNIDADES UNITS(Q1,M,Q2,DEGREE,Q3,DEGREE,Q4,DEGREE,Q5,DEGREE,Q6,DEGREE,Q7,DEGREE,U1,M/S,U 2,RAD/S,U3,RAD/S) UNITS(U4,RAD/S,U5,RAD/S,U6,RAD/S,U7,RAD/S,T,S) ! PROGRAMA EM FORTRAN CODE(CLIO63,SUBS) APÊNDICE 2 – Listagem do programa FORTRAN para a versão Clio 6.3 C THE NAME OF THIS PROGRAM IS CLIO63.FOR C CREATED BY AUTOLEV ON 10-18-2004 AT 08:48:50 C IMPLICIT DOUBLE PRECISION (A-Z) INTEGER JLOOP,NSTEPS,NCUTS,NEQS,ILOOP,COUNTER,NPSTEP LOGICAL STPSZ EXTERNAL EQNS CHARACTER MSG(75) DIMENSION U(14) COMMON/CZEES/Z(444) COMMON/CPAR/IS1,IS2,IS3,ID,IP,MS,MD,MP,M9,M11,PI,DEGTORAD,RADTODEG & ,DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L1,L2,L3,L4,K & D,KT,KP,CD,CT,V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FATOR3,FATOR4 COMMON/CONT/ZZ1,ZZ2,ZZ3,ZZ4,DYA9,DYA10,DYA11,DYA12,FA1,FA2,FA3,FA4 & ,TORQ COMMON/DFQLST/T,STEP,RELERR,ABSERR,NCUTS,NEQS,STPSZ C OPEN(UNIT=11,FILE='CLIO63.IN ',STATUS='UNKNOWN') OPEN(UNIT=12,FILE='CLIO63.OU1',STATUS='UNKNOWN') OPEN(UNIT=13,FILE='CLIO63.OU2',STATUS='UNKNOWN') OPEN(UNIT=14,FILE='CLIO63.OU3',STATUS='UNKNOWN') OPEN(UNIT=31,FILE='CLIO63.CO1',STATUS='UNKNOWN') OPEN(UNIT=32,FILE='CLIO63.CO2',STATUS='UNKNOWN') OPEN(UNIT=33,FILE='CLIO63.CO3',STATUS='UNKNOWN') PI = 4.0D0*DATAN(1.0D0) DEGTORAD = PI/180.0D0 RADTODEG = 1.0D0/DEGTORAD WRITE(*,6001) C ************************************************************************ * * * NOTE REGARDING INPUT AND OUTPUT DATA FILES * * * * The user must supply an input data file to this program. The * * file must be named FILENAME.IN , where FILENAME is obtained from * * the first line of this program. The data must be arranged in * * accordance with the READ statements that immediately follow this * * NOTE. * * * * The output from the program is sent to data files whose names * * appear on the screen at the completion of each run. The first * * column in each output data file contains the time T, * * running from zero to TMAX in increments of PSTEP. TMAX, PSTEP, * * and STEP are input from the terminal by the user at runtime, STEP * * being the initial integration stepsize, a number that should be * * chosen to be less than or equal to PSTEP. The terminal also * * prompts the user for a message identifying the run. This message * * is printed on each of the output files. Output files ending * * in .OUn contain time-histories of generalized speeds and * * generalized coordinates; files ending in .NRG contain kinetic * * energy, potential energy, and total energy time-histories; files * * ending in .H contain angular momentum time-histories; files ending * * in .COn contain time-histories of quantities appearing as * * arguments in CONTROLS commands; files ending in .SPn contain * * time-histories of SPECIFIED variables; and files ending in .AUn * * contain time-histories of force and/or torque measure numbers * * corresponding to AUXILIARY generalized speeds. * * * ************************************************************************ C READ(11,*) DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L1,L2 & ,L3,L4,KD,KT,KP,CD,CT,V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FATOR3 & ,FATOR4 READ(11,*) MS,MD,MP,M9,M11 READ(11,*) IS1,IS2,IS3 READ(11,*) ID READ(11,*) IP READ(11,*) U(1),U(2),U(3),U(4),U(5),U(6),U(7) READ(11,*) Q1,Q2,Q3,Q4,Q5,Q6,Q7 C C WRITE(* ,6002) READ( * ,6003) (MSG(ILOOP),ILOOP = 1,75) WRITE(* ,6009) READ(*,*) TMAX,PSTEP,STEP0 NPSTEP = IDINT((PSTEP-1.D-8)/STEP0 + 1) STEP = PSTEP/NPSTEP C WRITE(* ,6012) WRITE(* ,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(12,6101) WRITE(12,6010) (MSG(ILOOP),ILOOP = 1,75) 1 WRITE(13,6102) WRITE(13,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(14,6103) WRITE(14,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(31,6201) WRITE(31,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(32,6202) WRITE(32,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(33,6203) WRITE(33,6010) (MSG(ILOOP),ILOOP = 1,75) WRITE(* ,6011) WRITE(12,6011) C WRITE(* ,6500) DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L & 1,L2,L3,L4,KD,KT,KP,CD,CT,V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FA & TOR3,FATOR4 WRITE(12,6500) DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L & 1,L2,L3,L4,KD,KT,KP,CD,CT,V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FA & TOR3,FATOR4 WRITE(* ,6512) IS1,IS2,IS3 WRITE(12,6512) IS1,IS2,IS3 WRITE(* ,6514) ID WRITE(12,6514) ID WRITE(* ,6516) IP WRITE(12,6516) IP WRITE(* ,6600) MS,MD,MP,M9,M11 WRITE(12,6600) MS,MD,MP,M9,M11 WRITE(* ,6601) U(1),U(2),U(3),U(4),U(5),U(6),U(7) WRITE(12,6601) U(1),U(2),U(3),U(4),U(5),U(6),U(7) WRITE(* ,6602) Q1,Q2,Q3,Q4,Q5,Q6,Q7 WRITE(12,6602) Q1,Q2,Q3,Q4,Q5,Q6,Q7 WRITE(* ,6006) TMAX,PSTEP,STEP,STEP0 WRITE(12,6006) TMAX,PSTEP,STEP,STEP0 C U(8) = Q1 U(9) = Q2 U(10) = Q3 U(11) = Q4 U(12) = Q5 U(13) = Q6 U(14) = Q7 C WRITE(* ,6701) WRITE(12,6701) WRITE(13,6702) WRITE(14,6703) WRITE(31,6751) WRITE(32,6752) WRITE(33,6753) C NEQS = 14 NCUTS = 20 T = 0.0 RELERR = 1.0D-8 ABSERR = 1.0D-8 STPSZ = .FALSE. NSTEPS = IDINT(TMAX/STEP+0.1)+1 C COUNTER = 0 C C DO 1000 JLOOP = 1 , NSTEPS CALL ZEES(T,U) C IF (COUNTER.EQ.NPSTEP.OR.COUNTER.EQ.0) THEN WRITE(12,6005) T,U(1),U(2),U(3),U(4),U(5) WRITE(* ,6005) T,U(1),U(2),U(3),U(4),U(5) WRITE(13,6005) T,U(6),U(7),U(8),U(9),U(10) WRITE(14,6005) T,U(11),U(12),U(13),U(14) CALL CNTRL(T,U) WRITE(31,6005) T,ZZ1,ZZ2,ZZ3,ZZ4,DYA9 WRITE(32,6005) T,DYA10,DYA11,DYA12,FA1,FA2 WRITE(33,6005) T,FA3,FA4,TORQ COUNTER = 0 ENDIF C COUNTER = COUNTER + 1 IF (JLOOP.EQ.NSTEPS) GO TO 1000 CALL DEQS(EQNS,U,*99) 1000 CONTINUE C WRITE(*,6999) 2 C STOP 99 WRITE(*,6004) C 6001 FORMAT(/1X,'SYSTEM PARAMETERS AND INITIAL CONDITIONS'/ & 2X,'ARE NOW BEING READ FROM THE INPUT FILE'/) 6002 FORMAT(1X,'INPUT A DESCRIPTION OF THIS RUN'/) 6003 FORMAT(75A1) 6004 FORMAT(1X,'STEPSIZE HALVED TOO MANY TIMES'/) 6005 FORMAT(6(1X,1PE12.5)) 6006 FORMAT(11X,'TMAX = ',1PE12.5,' S'/10X,'PSTEP = ',1PE12.5,' S'/11X, &'STEP = ',1PE12.5,' S (USER INPUT VALUE = ',1PE12.5,' S)'//) 6007 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'HN1',10X,'HN2',10X,' &HN3',10X,'HN',/6X,'(S)',8X,'(UNITS)',6X,'(UNITS)',6X,'(UNITS)',6X, &'(UNITS)',/) 6008 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'KE',11X,'PE',9X,'KE &+ PE',/6X,'(S)',8X,'(UNITS)',6X,'(UNITS)',6X,'(UNITS)',/) 6009 FORMAT(/1X,'INPUT TMAX, PSTEP, STEP '/ &1X,'============================================='/ &1X,'| TMAX : FINAL TIME |'/ &1X,'| PSTEP: TIME INTERVAL FOR PRINTING |'/ &1X,'| STEP : MAXIMUM INTEGRATION TIME STEP |'/ &1X,'============================================='/) 6010 FORMAT(1X,'*** ',75A1) 6011 FORMAT(//1X,'SYSTEM PARAMETERS'/) 6012 FORMAT(1X,'OUTPUT FROM PROGRAM CLIO63.FOR'//) 6101 FORMAT(1X,'FILE: CLIO63.OU1 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6102 FORMAT(1X,'FILE: CLIO63.OU2 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6103 FORMAT(1X,'FILE: CLIO63.OU3 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6201 FORMAT(1X,'FILE: CLIO63.CO1 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6202 FORMAT(1X,'FILE: CLIO63.CO2 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6203 FORMAT(1X,'FILE: CLIO63.CO3 (OUTPUT FROM PROGRAM CLIO63.FOR)'//) 6500 FORMAT(12X,'DX1 = ',1PE12.5,' UNITS'/12X,'DY1 = ',1PE12.5,' UNITS' &/12X,'DZ1 = ',1PE12.5,' UNITS'/12X,'DX3 = ',1PE12.5,' UNITS'/12X,' &DY3 = ',1PE12.5,' UNITS'/12X,'DZ3 = ',1PE12.5,' UNITS'/12X,'DX5 = &',1PE12.5,' UNITS'/12X,'DY5 = ',1PE12.5,' UNITS'/12X,'DZ5 = ',1PE1 &2.5,' UNITS'/12X,'DX7 = ',1PE12.5,' UNITS'/12X,'DY7 = ',1PE12.5,' &UNITS'/12X,'DZ7 = ',1PE12.5,' UNITS'/14X,'G = ',1PE12.5,' UNITS'/1 &3X,'L1 = ',1PE12.5,' UNITS'/13X,'L2 = ',1PE12.5,' UNITS'/13X,'L3 = & ',1PE12.5,' UNITS'/13X,'L4 = ',1PE12.5,' UNITS'/13X,'KD = ',1PE12 &.5,' UNITS'/13X,'KT = ',1PE12.5,' UNITS'/13X,'KP = ',1PE12.5,' UNI &TS'/13X,'CD = ',1PE12.5,' UNITS'/13X,'CT = ',1PE12.5,' UNITS'/14X, &'V = ',1PE12.5,' UNITS'/14X,'L = ',1PE12.5,' UNITS'/13X,'X0 = ',1P &E12.5,' UNITS'/11X,'AMPL = ',1PE12.5,' UNITS'/9X,'LAMBDA = ',1PE12 &.5,' UNITS'/13X,'NL = ',1PE12.5,' UNITS'/9X,'FATOR1 = ',1PE12.5,' &UNITS'/9X,'FATOR2 = ',1PE12.5,' UNITS'/9X,'FATOR3 = ',1PE12.5,' UN &ITS'/9X,'FATOR4 = ',1PE12.5,' UNITS'/) 6512 FORMAT(12X,'IS1 = ',1PE12.5,' UNITS'/12X,'IS2 = ',1PE12.5,' UNITS' &/12X,'IS3 = ',1PE12.5,' UNITS'/) 6514 FORMAT(13X,'ID = ',1PE12.5,' UNITS'/) 6516 FORMAT(13X,'IP = ',1PE12.5,' UNITS'/) 6600 FORMAT(13X,'MS = ',1PE12.5,' UNITS'/13X,'MD = ',1PE12.5,' UNITS'/1 &3X,'MP = ',1PE12.5,' UNITS'/13X,'M9 = ',1PE12.5,' UNITS'/12X,'M11 &= ',1PE12.5,' UNITS'//) 6601 FORMAT(/1X,'INITIAL CONDITIONS'//10X,'U1(0) = ',1PE12.5,' M/S'/10X &,'U2(0) = ',1PE12.5,' RAD/S'/10X,'U3(0) = ',1PE12.5,' RAD/S'/10X,' &U4(0) = ',1PE12.5,' RAD/S'/10X,'U5(0) = ',1PE12.5,' RAD/S'/10X,'U6 &(0) = ',1PE12.5,' RAD/S'/10X,'U7(0) = ',1PE12.5,' RAD/S'/) 6602 FORMAT(10X,'Q1(0) = ',1PE12.5,' M'/10X,'Q2(0) = ',1PE12.5,' DEGREE &'/10X,'Q3(0) = ',1PE12.5,' DEGREE'/10X,'Q4(0) = ',1PE12.5,' DEGREE &'/10X,'Q5(0) = ',1PE12.5,' DEGREE'/10X,'Q6(0) = ',1PE12.5,' DEGREE &'/10X,'Q7(0) = ',1PE12.5,' DEGREE'//) 6701 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'U1',11X,'U2',11X,'U3 &',11X,'U4',11X,'U5',/6X,'(S)',9X,'(M/S)',7X,'(RAD/S)',6X,'(RAD/S)' &,6X,'(RAD/S)',6X,'(RAD/S)',/) 6702 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'U6',11X,'U7',11X,'Q1 &',11X,'Q2',11X,'Q3',/6X,'(S)',8X,'(RAD/S)',6X,'(RAD/S)',8X,'(M)',7 &X,'(DEGREE)',5X,'(DEGREE)',/) 6703 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'Q4',11X,'Q5',11X,'Q6 &',11X,'Q7',/6X,'(S)',7X,'(DEGREE)',5X,'(DEGREE)',5X,'(DEGREE)',5X, &'(DEGREE)',/) 6751 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'ZZ1',10X,'ZZ2',10X,' &ZZ3',10X,'ZZ4',9X,'DYA9',/6X,'(S)',8X,'(UNITS)',6X,'(UNITS)',6X,'( &UNITS)',6X,'(UNITS)',6X,'(UNITS)',/) 6752 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',10X,'DYA10',8X,'DYA11',8X &,'DYA12',9X,'FA1',10X,'FA2',/6X,'(S)',8X,'(UNITS)',6X,'(UNITS)',6X &,'(UNITS)',6X,'(UNITS)',6X,'(UNITS)',/) 6753 FORMAT(//1X,'SIMULATION RESULTS'//7X,'T',11X,'FA3',10X,'FA4',9X,'T &ORQ',/6X,'(S)',8X,'(UNITS)',6X,'(UNITS)',6X,'(UNITS)',/) 6999 FORMAT(//1X,'OUTPUT IS ON FILES: ','CLIO63.OU1'/22X,'CLIO63.OU2'/ &22X,'CLIO63.OU3'/22X,'CLIO63.CO1'/22X,'CLIO63.CO2'/22X,'CLIO63.CO3 &'/22X,/) 3 END C C C SUBROUTINE EQNS(T,U,UDOT) IMPLICIT DOUBLE PRECISION (A-Z) DIMENSION U(14),UDOT(14),COEF(7,7),RHS(7) COMMON/CZEES/Z(444) COMMON/CPAR/IS1,IS2,IS3,ID,IP,MS,MD,MP,M9,M11,PI,DEGTORAD,RADTODEG & ,DX1,DY1,DZ1,DX3,DY3,DZ3,DX5,DY5,DZ5,DX7,DY7,DZ7,G,L1,L2,L3,L4,K & D,KT,KP,CD,CT,V,L,X0,AMPL,LAMBDA,NL,FATOR1,FATOR2,FATOR3,FATOR4 COMMON/CONT/ZZ1,ZZ2,ZZ3,ZZ4,DYA9,DYA10,DYA11,DYA12,FA1,FA2,FA3,FA4 & ,TORQ C C CALL ZEES(T,U) C Q1 = U(8) Q2 = U(9) Q3 = U(10) Q4 = U(11) Q5 = U(12) Q6 = U(13) Q7 = U(14) C C CALL CNTRL(T,U) C C S1 = DSIN(Q2) C1 = DCOS(Q2) C COEF(1,1) = -Z(413) COEF(1,2) = Z(414) COEF(1,3) = Z(415) COEF(1,4) = Z(416) COEF(1,5) = -Z(417) COEF(1,6) = -Z(418) COEF(1,7) = -Z(419) COEF(2,1) = Z(414) COEF(2,2) = Z(421) COEF(2,3) = Z(422) COEF(2,4) = Z(423) COEF(2,5) = Z(424) COEF(2,6) = Z(425) COEF(2,7) = Z(426) COEF(3,1) = Z(415) COEF(3,2) = Z(428) COEF(3,3) = Z(429) COEF(3,4) = Z(430) COEF(3,5) = Z(431) COEF(3,6) = Z(432) COEF(3,7) = Z(433) COEF(4,1) = Z(416) COEF(4,2) = Z(423) COEF(4,3) = Z(430) COEF(4,4) = Z(435) COEF(4,5) = 0.0 COEF(4,6) = 0.0 COEF(4,7) = 0.0 COEF(5,1) = -Z(417) COEF(5,2) = Z(424) COEF(5,3) = Z(431) COEF(5,4) = 0.0 COEF(5,5) = Z(437) COEF(5,6) = 0.0 COEF(5,7) = 0.0 COEF(6,1) = -Z(418) COEF(6,2) = Z(439) COEF(6,3) = Z(432) COEF(6,4) = 0.0 COEF(6,5) = 0.0 COEF(6,6) = Z(440) COEF(6,7) = 0.0 COEF(7,1) = -Z(419) COEF(7,2) = Z(442) COEF(7,3) = Z(433) COEF(7,4) = 0.0 COEF(7,5) = 0.0 COEF(7,6) = 0.0 COEF(7,7) = Z(443) C 4 RHS(1) = ((-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L1*Z(19)+Q1)-ZZ2)*KP+((-D & X5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)-L1*Z(9)+Q1)-ZZ1)*KP+((DX7*S1+DY7*Z(3 & )-DZ7+DZ7*Z(4)+L4*Z(41)+Q1)-ZZ4)*KP+((DX7*S1-DY7*Z(3)-DZ7+DZ7*Z( & 4)+L4*Z(29)+Q1)-ZZ3)*KP-2*G*M11-2*G*M9-2*G*MD-2*G*MP-G*MS-Z(420) RHS(2) = -((C1*DX3*U(2)+DY3*Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)+D & Y7*Z(187)-DZ7*Z(188)+L3*Z(248)+U(1)))*CT*Z(334)-((C1*DX3*U(2)-DY & 3*Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)-DY7*Z(187)-DZ7*Z(188)+L3* & Z(229)+U(1)))*CT*Z(336)-(-((-C1*DX1*U(2)+DY1*Z(187)+DZ1*Z(188)+U & (1))-(-C1*DX5*U(2)+DY5*Z(187)-DZ5*Z(188)+L2*Z(210)+U(1)))*CD-((- & DX1*S1+DY1*Z(3)+DZ1-DZ1*Z(4)+Q1)-(-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+ & L2*Z(19)+Q1))*KD)*Z(338)+(-((-C1*DX1*U(2)-DY1*Z(187)+DZ1*Z(188)+ & U(1))-(-C1*DX5*U(2)-DY5*Z(187)-DZ5*Z(188)-L2*Z(191)+U(1)))*CD-(( & -DX1*S1-DY1*Z(3)+DZ1-DZ1*Z(4)+Q1)-(-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4) & -L2*Z(9)+Q1))*KD)*Z(340)+(-((-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L1*Z( & 19)+Q1)-ZZ2)*KP+G*M9)*Z(125)-(-((-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)-L & 1*Z(9)+Q1)-ZZ1)*KP+G*M9)*Z(100)-(-((DX7*S1+DY7*Z(3)-DZ7+DZ7*Z(4) & +L4*Z(41)+Q1)-ZZ4)*KP+G*M11)*Z(174)-(-((DX7*S1-DY7*Z(3)-DZ7+DZ7* & Z(4)+L4*Z(29)+Q1)-ZZ3)*KP+G*M11)*Z(149)-G*MD*Z(109)+G*MD*Z(133)- & G*MP*Z(158)-G*MP*Z(182)+Z(427) RHS(3) = ((C1*DX3*U(2)+DY3*Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)+DY & 7*Z(187)-DZ7*Z(188)+L3*Z(248)+U(1)))*CT*Z(343)+((C1*DX3*U(2)-DY3 & *Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)-DY7*Z(187)-DZ7*Z(188)+L3*Z & (229)+U(1)))*CT*Z(346)+(-((-C1*DX1*U(2)+DY1*Z(187)+DZ1*Z(188)+U( & 1))-(-C1*DX5*U(2)+DY5*Z(187)-DZ5*Z(188)+L2*Z(210)+U(1)))*CD-((-D & X1*S1+DY1*Z(3)+DZ1-DZ1*Z(4)+Q1)-(-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L & 2*Z(19)+Q1))*KD)*Z(349)-(-((-C1*DX1*U(2)-DY1*Z(187)+DZ1*Z(188)+U & (1))-(-C1*DX5*U(2)-DY5*Z(187)-DZ5*Z(188)-L2*Z(191)+U(1)))*CD-((- & DX1*S1-DY1*Z(3)+DZ1-DZ1*Z(4)+Q1)-(-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)- & L2*Z(9)+Q1))*KD)*Z(352)-(-((-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L1*Z(1 & 9)+Q1)-ZZ2)*KP+G*M9)*Z(126)+(-((-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)-L1 & *Z(9)+Q1)-ZZ1)*KP+G*M9)*Z(101)+(-((DX7*S1+DY7*Z(3)-DZ7+DZ7*Z(4)+ & L4*Z(41)+Q1)-ZZ4)*KP+G*M11)*Z(175)+(-((DX7*S1-DY7*Z(3)-DZ7+DZ7*Z & (4)+L4*Z(29)+Q1)-ZZ3)*KP+G*M11)*Z(150)+G*MD*Z(110)-G*MD*Z(134)+G & *MP*Z(159)+G*MP*Z(183)+Z(434) RHS(4) = -(-((-C1*DX1*U(2)-DY1*Z(187)+DZ1*Z(188)+U(1))-(-C1*DX5*U( & 2)-DY5*Z(187)-DZ5*Z(188)-L2*Z(191)+U(1)))*CD-((-DX1*S1-DY1*Z(3)+ & DZ1-DZ1*Z(4)+Q1)-(-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)-L2*Z(9)+Q1))*KD) & *Z(94)+(-((-DX5*S1-DY5*Z(3)-DZ5+DZ5*Z(4)-L1*Z(9)+Q1)-ZZ1)*KP+G*M & 9)*Z(102)+G*MD*Z(111)+Z(436) RHS(5) = (-((-C1*DX1*U(2)+DY1*Z(187)+DZ1*Z(188)+U(1))-(-C1*DX5*U(2 & )+DY5*Z(187)-DZ5*Z(188)+L2*Z(210)+U(1)))*CD-((-DX1*S1+DY1*Z(3)+D & Z1-DZ1*Z(4)+Q1)-(-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L2*Z(19)+Q1))*KD) & *Z(119)-(-((-DX5*S1+DY5*Z(3)-DZ5+DZ5*Z(4)+L1*Z(19)+Q1)-ZZ2)*KP+G & *M9)*Z(127)-G*MD*Z(135)+Z(438) RHS(6) = -((C1*DX3*U(2)-DY3*Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)-D & Y7*Z(187)-DZ7*Z(188)+L3*Z(229)+U(1)))*CT*Z(143)-(-((DX7*S1-DY7*Z & (3)-DZ7+DZ7*Z(4)+L4*Z(29)+Q1)-ZZ3)*KP+G*M11)*Z(151)-G*MP*Z(160)+ & KT*Q6+Z(441) RHS(7) = -((C1*DX3*U(2)+DY3*Z(187)-DZ3*Z(188)+U(1))-(C1*DX7*U(2)+D & Y7*Z(187)-DZ7*Z(188)+L3*Z(248)+U(1)))*CT*Z(168)-(-((DX7*S1+DY7*Z & (3)-DZ7+DZ7*Z(4)+L4*Z(41)+Q1)-ZZ4)*KP+G*M11)*Z(176)-G*MP*Z(184)+ & KT*Q7+Z(444) C CALL UNCUPL(7,COEF,RHS,UDOT) C C C U8 IS DEFINED TO BE Q1 UDOT(8) = U(1) C C U9 I