edit · print · PDF

Please note that all the SIEpedia's articles address specific issues or questions raised by IAC users, so they do not attempt to be rigorous or exhaustive, and may or may not be useful or applicable in different or more general contexts.

La necesidad de hacer Fortran 95 compatible con las versiones anteriores hace que haya maneras muy diferentes de hacer cada cosa. Consecuencia de esto es la facilidad de caer en malas prácticas que hagan nuestro código poco portable, o incluso inservible ante un simple cambio de compilador o de cpu. Que lo hayamos probado con éxito en nuestra máquina no garantiza que funcione en otra.

Esta sección pretende mostrar "recetas" de trabajo que permitan programar de modo sólido y eviten las trampas que, en aras de la compatibilidad con estándares antiguos, se esconden en el estándar Fortran95.

Declaración y asignación de valores de coma flotante

Si usamos las viejas fórmulas para declarar simple y doble precisión ("REAL v" y "DOUBLE PRECISION v") después será casi imposible cambiar la precisión en todo el código sin tener problemas: para pasar de simple a doble, no bastará con cambiar la declaración, sino que tendremos que cambiar UNA A UNA todas las asignaciones del código. Un ejemplo de la forma correcta es este:

Primero asignamos al entero SS (por ejemplo) el valor del parámetro kind que el sistema asigna para una mantisa de (por ejemplo) seis cifras:

  INTEGER, PARAMETER :: SS = SELECTED_REAL_KIND(6)

Los flotantes los declaramos entonces así:

  REAL(KIND = SS) :: v, w

Los valores numéricos en asignaciones dentro del código, precisan entonces un sufijo _SS así:

  v = 5.050_SS
  w = 1.0E-3_SS
  (...)

Esta es la mejor manera de asegurar que, en cualquier sistema, tendremos una mantisa de al menos 6 dígitos (mas otro oculto, pero eso ya es otro tema) y de que, si queremos cambiar la precisión en todo el código, lo haremos de golpe, solo con cambiar el valor X en SELECTED_REAL_KIND(x)

Importante: No hay que olvidar el sufijo (_SS en el ejemplo) en las asignaciones; es un error muy común.

Podemos tener varios tipos de precisión, así como declarar enteros más cortos (para ganar velocidad). Ejemplito:

 PROGRAM ucho
  IMPLICIT NONE
  !
  INTEGER, PARAMETER :: SS = SELECTED_REAL_KIND(6)     ! mantisa de al menos 6 cifras
  INTEGER, PARAMETER :: DD = SELECTED_REAL_KIND(15)    ! mantisa de al menos 15 cifras
  INTEGER, PARAMETER :: QQ = SELECTED_REAL_KIND(33)    ! mantisa de al menos 33 cifras
  INTEGER, PARAMETER :: II = SELECTED_INT_KIND(3)      ! enteros de tres dígitos, más cortos pero más rápidos
  !
  REAL(KIND = ss)     :: s
  REAL(KIND = DD)     :: d
  REAL(KIND = QQ)     :: q
  INTEGER (KIND = II) :: i
  !
  s = ACOS(-1.0_SS)
  d = ACOS(-1.0_DD)
  q = ACOS(-1.0_QQ)
  !
  PRINT*, s     ! 3.141593
  PRINT*, d     ! 3.14159265358979
  PRINT*, q     ! 3.1415926535897932384626433832795028
  !
  END PROGRAM ucho

No caigais tampoco en cosas como "REAL(KIND = 8)", porque ese 8 puede tener un significado diferente en otra máquina. Sólo es seguro obtener primero el valor mediante la función SELECTED_REAL_KIND(x)

Hay que huir también de formas fuera del estándar, aunque aún sean muy populares, como "REAL*8", etc.

Operaciones lógicas con valores de coma flotante

Los números con decimales que la máquina nos presenta, son una deferencia que tiene hacia los humanos, incapaces de leer la verdadera cantidad con la que la máquina está trabajando: una cadena binaria que representa una fracción entera.

Esa fracción entera es en general solo una aproximación, dentro de las N+1 cifras significativas del modelo declarado (N=6 en simple precisión) pero no más allá de la cifra N+1. Por ejemplo, el número 0.01 en realidad se almacena como una fracción que "en verdad" es igual a 0.009999999776482...

¿Qué nos importa eso? Nos importa, porque la máquina solo ve la fracción entera. Por eso, las operaciones lógicas protagonizadas por comparaciones directas entre reales, son una práctica suicida, porque no sabemos realmente cuál es el valor que está detras. He aquí una recetilla para sortear esto:

Como comparar dos flotantes
Hay que comparar su diferencia con un valor muy pequeño. A su vez, hay que definir ese valor como "pequeño" en relación a las propias cantidades a comparar, lo que dependerá del tipo de flotante que hayamos declarado; para ello tenemos la función intrínseca EPSILON(x). Un ejemplito:

  integer, parameter :: p = selected_real_kind(6)
  real(kind=p) :: a, b
  ! (...)
  if ( abs(a-b) <= abs(a+b)*epsilon(a) ) then
    print*, "los numericos son iguales"
  else if ( (a-b) > abs(a+b)*epsilon(a) ) then
    print*, "A es más grande que B"
  else if ( (b-a) > abs(a+b)*epsilon(a) ) then
    print*, "B es más grande que A"
  endif

Para quien lo quiera, aquí está definido el procedimiento en una función, que compara dos flotantes y devuelve -1, 1 ó 0:

  function compara(x1,x2)
    !! Devuelve: 0 si x1 y x2 son iguales,
    !!           1 si x1 > x2
    !!          -1 si x1 < x2
    !! Requiere haber declarado en el módulo el kind de los flotantes en p, por ejemplo:
    !! integer, parameter, public :: p = selected_real_kind(6)
    !! (si no la usamos en un módulo, declarar p sin el atributo PUBLIC)
    real(kind=p), intent(in)  :: x1, x2
    integer :: compara
    if ( abs(x1-x2) <= abs(x1+x2)*epsilon(x1) ) then
      compara = 0
    else if ( (x1-x2) > abs(x1+x2)*epsilon(x1) ) then
      compara = 1
    else if ( (x2-x1) > abs(x1+x2)*epsilon(x1) ) then
      compara = -1
    endif
  end function compara

Por ejemplo,

 print*, compara(3.001, 3.002)

saca por pantalla el entero -1

También hay que tener mucho cuidado al truncar resultados. Aunque "matemáticamente" algo deba dar 1.000, si "por dentro" es 0.9999... y lo truncamos directamente con la función int(), nos quedará cero.

Escribir correctamente a disco/ pantalla los valores de coma flotante

No os comáis más la cabeza con esto: estos campos de formato funcionan muy bien para escribir valores de coma flotante:

  es13.6e2
  es23.15e3

Para simple y doble precisión, respectivamente, esto es, kind = selected_real_kind(6) y kind = selected_real_kind(15)

Centrémonos en el primero de ellos para aclarar el por qué: la mantisa de simple precisión, de parámetro kind obtenido mediante selected_real_kind(6), tiene al menos seis cifras significativas, pero en realidad trabaja con dos decimales más, uno de ellos exacto y el otro aproximado. Podemos contar entonces con un máximo de seis decimales (además de la cifra de los enteros). Por otro lado, el exponente tendrá un máximo de dos dígintos. Junto con el punto decimal, el posible signo y demás, el campo mínimo ha de tener 13 huecos.

No tiene sentido sacar más decimales. Tampoco es buena práctica usar los formatos tipo f8.5 y demás. La única manera de escribir la mantisa completa y sin sacar decimales espúreos es usando esos formatos.

Otra cosa es querer sacar mantisas más cortas, para redondear los resultados. Por ejemplo, es12.5e2 con un decimal menos.

¿Cómo me aseguro de que mi estilo de programación no cae en prácticas problemáticas y/o obsoletas?

En primera instancia, ayuda bastante el pasarle al compilador opciones para que trate de avisarnos sobre las malas prácticas, y para que nos prohíba usar extensiones al estándar. A mí me funciona bien esto:

Para ifort de Intel:

  ifort -e95 -warn -error-limit 1

Para g95:

  g95 -std=f95 -pedantic -Wall -Wextra -fone-error

... donde además les decimos que nos den un solo mensaje de error cada vez.

Aún así, nada nos dice cuál es la manera óptima de hacer cada cosa, y es un poco agotador andar siempre tomando decisiones, tratando de averiguar cuál es la mejor forma para algo: ¿Subrutinas en módulos o adjuntas con interface explícita? ¿tabular dos espacios con TR2 ó 2x? ¿Por qué hay especificadores redundantes en muchas sentencias, como END, ERR e IOSTAT en la sentencia OPEN, y cuál es el bueno y cuál el que se mantiene por compatibilidad pero es peor? Para evitar esos dolores de cabeza e ir directamente a lo óptimo, el compilador g95 nos ofrece a los puristas la opción destroyer:

  -std=F

Le estamos diciendo que se restrinja a un subconjunto dentro del estándar Fortran 95, llamado F. Unos tipos en la Royal University of Nosedonde, preocupados por el problema de la ambigüedad en Fortran95, seleccionaron un subconjunto de reglas de sintaxis dentro del estándar Fortran 95, que dejaran de lado todas las prácticas peligrosas y anticuadas.

Acostumbrarse a F es un poco agobiante al principio, pero te hace ganar solidez y dominio de las capacidades de Fortran 95 muy rápidamente. Como han seleccionado una y solo una forma de hacer cada cosa, en pocas semanas llegas al dominio total de Fortran 95.

La sintaxis de F se especifica en: http://www.fortran.com/F/F_bnf.html

Además, se use F o no, estaría muy bien echar un ojo a este documento, con recomendaciones de buenas prácticas: "Standards, Guidelines and Recommendations for Writing Fortran 95 Code" http://projects.osd.noaa.gov/spsrb/standards_docs/Fortran95_standard_rev22Jun2009.pdf

Ayuda online para la sintaxis, y buena documentación

Es muy útil tener este link siempre a mano, para consultas rápidas:

"Compact Fortran 95 Language Summary" http://www.cs.umbc.edu/~squire/fortranclass/summary.shtml

Y, para consultas más serias, lo mejor es dejarse de rodeos y acudir directamente al documento del estándar ANSI. Porque en los libros, los autores suelen seleccionar unas cosas y callar otras, para resumir. Solo hay un sitio con toda la información exhaustiva, aunque su formato sea un poco árido (parece el código penal...) y es el estándar. Dos links:

http://j3-fortran.org/doc/standing/archive/007/97-007r2/pdf/97-007r2.pdf
ftp://ftp.nag.co.uk/sc22wg5/N1151-N1200/N1191.pdf



Eduardo Guerras (eguerras@iac.es)

edit · print · PDF
Page last modified on August 06, 2010, at 02:35 PM