Cálculo Numérico com Fortran – Parte 1

Depois de um longo jejum de escritos no meu blog, resolvi começar uma nova empreitada, a intenção é aprender Fortran aplicado ao cálculo numérico.

Para os estudos de caso, seguirei a tradução do livro Numerical Methods and FORTRAN Programming (Wiley 64). Porem o livro foi escrito em Fortran IV então fiz um upgrade no código para a versão 90, e como sempre não me limitarei a copiar o código e posso terminar fazendo alguma implementação própria.

Para quem não sabe o que é Cálculo Numérico http://www.alunos.eel.usp.br/numerico/notasDeAula/intro.pdf

Para quem não sabe o que é Fortran http://pt.wikipedia.org/wiki/Fortran

Aqui segue um tutorial para instalação do gfortran e gnuplot http://www.fsc.ufsc.br/~esalcedo/FisComp/Notas/InstalarNowindows/Fortran/f90windows.html

Estou estudando o livro na medida que escrevo por isso, os primeiros algoritmos podem ser bem toscos, pois não domino a linguagem (AINDA) :p

O primeiro estudo de caso proposto por Wiley é o calculo das raízes de um polinômio de 4º grau pelo método da bisseção.

O metodo da bissecção consiste em procurar intervalos valores para x0 e x1, que tornem verdadeiras a afirmativa abaixo F(x0) * F(x1)<=0 , e reduzir esse intervalo cada vez mais até que o valor de F(x) seja o mais próximo ou igual a 0

Para entender melhor usaremos a função y=x^4-9x^3-2x^2+120x-130

Vejamos o gráfico dessa função:

No intervalo de -10 a 10 temos 4 raízes* visualmente podemos distinguir que a primeira esta entre -4 e -3 a segunda entre 1 e 2 a terceira entre 3 e 4 e a quarta entre 7 e 8. A bissecção consiste em pegar esses intervalos e reduzi-los, tomamos como exemplo o intervalo entre -4 e -3, pegamos então o ponto médio entre eles que chamamos de x,  calculamos então o valor de F(x), existem 3 possibilidades, ou F(x) = 0 caso que x é uma raiz, ou F(x)>0 caso em que o a raiz esta entre x e -3, ou F(x) < 0 caso que a raiz esta entre -4 e x.

* Raízes de uma função são valores de X que tornam F(x) = 0

O fluxograma abaixo, demonstra a bissecção dessa função:

Implementando em Fortran:

program teste

REAL*4 XL/0/
 REAL*4 XR/0/
 REAL*4 YR/0/
 REAL*4 YL/0/
 REAL*4 XS/0/
 REAL*4 X/0/
 REAL*4 Y/0/

 WRITE (*,*) "Estudo de Caso 1: Determinacao de raizes por bisseccao"
 PRINT *,"Entre um Byte e Outro - Aprendendo Fortran e Calculo Numerico"
 PRINT *,"Baseado na Traducao do livro Numerical Methods and Fortran Programming (Wiley, 1964)"
 PRINT *,"Adaptação do codigo para Fortran 90 por Anthony Collucci"
 PRINT *,"Esse algoritmo supoe que todas as raizes da funcao estão no intervalo entre -10 e 10"
 Print *," e que entre cada raiz a um intervalo nao inferior a 1"

 !Declarando o Intervalo de pesquisa o XL é o limite inferior.
 XL = - 10
 XR = XL + 1

 PRINT *," Xl ",","," Xr "
 DO
 ! Calcula o valor da função para XL E XR
 Yl = F(XL)
 YR = F(XR)

 PRINT *,XL,",",XR

 ! Caso As raizes estejam entre XL e XR, o produto de YR, e YL sera negativo pois o gráfico da função
 ! cruza o eixo X.

 IF (YR*YL.LT.0) THEN
 !SALVA O VALOR DE XR
 XS = XR

 DO

 !Define o ponto x entre XR e XL
 X = (XR+XL)/2
 Y = F (X)

 IF (YL*Y.EQ.0) THEN
 !Se o produto e 0 então a raiz e x, pois senao haveriamos encontrado ela no if anterior ao loop
 PRINT *,"RAIZ E", X
 EXIT
 ELSEIF (YL*Y.LT.0) THEN
 !Se for menor que 0 significa que a raiz esta entre x e xl
 XR = X
 YR = Y
 ELSEIF (YL*Y.GT.0) THEN
 !Se for maior que 0 significa que esta entre x e xr
 XL = X
 YL = Y
 ENDIF

 ! Se a precisao desejada foi alcancada entao paramos
 IF (XR-XL.LT.0.000005) THEN
 PRINT *,"RAIZ ESTA ENTRE",XL, ",", XR
 !Restaura o valor de XR
 XR = XS
 EXIT
 ENDIF
 END DO
 ENDIF

 ! Aqui definimos que procuraremos uma raiz até o XR = 10.
 IF (XR.GT.10) EXIT

 ! Incrementamos o valor de XR e XL.
 XL = XR
 XR = XL+1

 END DO

 STOP

 CONTAINS

 FUNCTION F(X)
 ! Função que criamos para calcular o valor da função matematica proposta.
 REAL*4, INTENT (IN) :: X
 REAL*4 Y/0/
 Y = X**4 - 9*X**3 - 2*X**2 +120*X - 130
 F=Y
 RETURN
 END FUNCTION F

 end program teste

Veja como fica a saída do programa:

 

 

Links Uteis (Contendo as referencias para a sintaxe):

http://www.fsc.ufsc.br/~esalcedo/FisComp/Notas/InstalarNowindows/Fortran/f90windows.html
http://nf.nci.org.au/training/FortranAdvanced/slides/slides.015.html
http://www.cs.mtu.edu/~shene/COURSES/cs201/NOTES/F90-Basics.pdf
http://www.inf.ufes.br/~thomas/fortran/tutorials/helder/fortran.pdf
http://www.dsc.ufcg.edu.br/~icc/Periodo-2011.2/matadic/Programacao%20em%20Linguagem%20FORTRAN.htm#P050
http://www.ews.uiuc.edu/~mrgates2/docs/fortran.html

Uma ideia sobre “Cálculo Numérico com Fortran – Parte 1”

Deixe uma resposta

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *

Please type the characters of this captcha image in the input box

Por favor, digite os caracteres desta imagem na caixa de entrada