Fortran [PDF]

Zitiervorschau

Dan RacoŃi Elemente de Fortran 95 Cuprins 0. Introducere 1.0 Constante şi variabile 1.1 OperaŃii algebrice simple 2.1 InstrucŃiunea IF 2.1.1 InstrucŃiunea IF aritmetic 2.1.2 InstrucŃiunea IF logic 2.1.3 InstrucŃiunea BLOCK IF 2.1.4 Compararea şirurilor de caractere 2.1.5 InstrucŃiunea select case 3. InstrucŃiunea DO 3.1 Forma 1 3.2 Forma 2 3.3 Forma 3 3.4 InstrucŃiunea DO WHILE 3.5 IteraŃii 3.5.1 Calculul radicalului 3.5.2 Rezolvarea ecuaŃiilor algebrice neliniare cu metoda Newton 3.5.3 Rezolvarea sistemelor algebrice neliniare de două ecuaŃii cu metoda Newton 3.5.4 Rezolvarea ecuaŃiilor algebrice polinomiale cu coeficienŃi complecşi cu metoda Newton 4. Dezvoltări în serie 5. Tablouri 5.1 Vectori 5.2 Matrice 5.3 Gauss Seidel 5.4 Alocarea dinamică de memorie (în execuŃie) 6. FuncŃii şi subprograme 6.1 Clauza contains. FuncŃii interne 6.2 FuncŃii externe 6.3 FuncŃii în complex 6.4 InstrucŃiunea external 6.5 Subrutine 6.5.1 Transmiterea tablourilor la subprograme 6.5.2 Rezolvarea sistemelor de ecuaŃii algebrice lineare 6.5.3 Calculul valorilor proprii 7. InstrucŃiunea COMMON 7.1 InstrucŃiunea COMMON blank 7.2 InstrucŃiunea COMMON etichetat 7.3 InstrucŃiunea BLOCK DATA 8. InstrucŃiunea INTERFACE 9. InstrucŃiunea MODULE 10. AplicaŃii 10.1 Calculul integralelor definite, QUADPACK 10.2 Integrarea ecuaŃiilor diferenŃiale ordinare, problema Cauchy 10.3 EcuaŃia Burgers omogenă

1

Dan RacoŃi Elemente de Fortran 95 10.4 EcuaŃia Poisson 10.5 EcuaŃia de potenŃial pe cerc 10.6 Rezolvarea sistemelor de ecuaŃii neliniare 0. Introducere Fortran este un limbaj de programare potrivit în special pentru calculule numerice şi calcule ştiinŃifice. Limbajul dezvoltat plecând din 1950 a fost utilizat extensiv în meteorologie, analiza structurilor cu metoda elementelor finite, computational fluid dynamics (CFD), computational physics and computational chemistry. Este limbajul de programare utilizat pe supercalculatoare. Numele limbajului provine de la FORmula TRANslating System. Versiunile succesive au adăugat procesarea şirurilor de caractere, block if, do enddo, do while (FOTRAN 77), programarea modulară, array sections, object-based programming (Fortran 90/95) şi object-oriented and generic programming (Fortran 2003). Prin evoluŃie limbajul a devenit deosebit de complex. În expunerea succintă care urmează s-au avut în vedere în special aspectele legate de metodele numerice şi CFD. Pentru compatibilitate cu bibliotecile de programe existente, NETLIB.ORG, Fortran 90 Codes, etc., se prezintă toată gama de instrucŃiuni incluzând instrucŃiunea COMMON şi BLOCK DATA. Programarea modulară (MODULE) permite scrierea de programe compacte şi bine structurate şi de aceea este recomandată ca tehnică de programare. Sunt prezentate o serie de exemple simple (tehnici numerice de bază) dar şi aplicaŃii mai complexe (ecuaŃii diferenŃiale cu derivate parŃiale 2D). Unele programe au fost păstrate în forma originală deoarece sunt algoritmi omologaŃi. Forma fixă a limbajului *.f, *.for are 80 de coloane cu un caracter de continuare în coloana 6 (diferit de zero) . Zona 1-5 este zonă etichetă. Zona 7-72 este zonă instrucŃiune. Comentariile încep cu litera C din coloana 1 sau cu !. Următorul program face parte din biblioteca pppack din NETIB. c ivex.f chapter iv. runge example, with cubic hermite interpolation c from * a practical guide to splines * by c. de boor integer i,istep,j,n,nm1 real aloger,algerp,c(4,20),decay,divdf1,divdf3,dtau,dx,errmax,g,h * ,pnatx,step,tau(20) data step, istep /20., 20/ g(x) = 1./(1.+(5.*x)**2) print 600 600 format(28h n max.error decay exp.//) decay = 0. do 40 n=2,20,2 c choose interpolation points tau(1), ..., tau(n) , equally c spaced in (-1,1), and set c(1,i) = g(tau(i)), c(2,i) = c gprime(tau(i)) = -50.*tau(i)*g(tau(i))**2, i=1,...,n. nm1 = n-1 h = 2./float(nm1) do 10 i=1,n tau(i) = float(i-1)*h - 1. c(1,i) = g(tau(i)) 10 c(2,i) = -50.*tau(i)*c(1,i)**2 c calculate the coefficients of the polynomial pieces c do 20 i=1,nm1

2

Dan RacoŃi Elemente de Fortran 95

20 c c

dtau = divdf1 divdf3 c(3,i) c(4,i)

tau(i+1) - tau(i) = (c(1,i+1) - c(1,i))/dtau = c(2,i) + c(2,i+1) - 2.*divdf1 = (divdf1 - c(2,i) - divdf3)/dtau = (divdf3/dtau)/dtau

estimate max.interpolation error on (-1,1). errmax = 0. do 30 i=2,n dx = (tau(i)-tau(i-1))/step do 30 j=1,istep h = float(j)*dx evaluate (i-1)st cubic piece

c c

pnatx = c(1,i-1)+h*(c(2,i-1)+h*(c(3,i-1)+h*c(4,i-1))) c 30

errmax = amax1(errmax,abs(g(tau(i-1)+h)-pnatx)) aloger = alog(errmax) if (n .gt. 2) decay = * (aloger - algerp)/alog(float(n)/float(n-2)) algerp = aloger 40 print 640,n,errmax,decay 640 format(i3,e12.4,f11.2) stop end

Forma free a limbajului are extensia *.f90. Caracterul de continuare este &, la sfârşitul liniei. Comentariile încep cu semnul de exclamare !. ! ivex.f90 ! chapter iv. runge example, with cubic hermite interpolation ! from * a practical guide to splines * by c. de boor integer i,istep,j,n,nm1 real aloger,algerp,c(4,20),decay,divdf1,divdf3,dtau,dx,errmax,g,h & ,pnatx,step,tau(20) data step, istep /20., 20/ g(x) = 1./(1.+(5.*x)**2) print 600 600 format(28h n max.error decay exp.//) decay = 0. do 40 n=2,20,2 ! choose interpolation points tau(1), ..., tau(n) , equally ! spaced in (-1,1), and set c(1,i) = g(tau(i)), c(2,i) = ! gprime(tau(i)) = -50.*tau(i)*g(tau(i))**2, i=1,...,n. nm1 = n-1 h = 2./float(nm1) do 10 i=1,n tau(i) = float(i-1)*h - 1. c(1,i) = g(tau(i)) 10 c(2,i) = -50.*tau(i)*c(1,i)**2 ! calculate the coefficients of the polynomial pieces ! do 20 i=1,nm1 dtau = tau(i+1) - tau(i) divdf1 = (c(1,i+1) - c(1,i))/dtau divdf3 = c(2,i) + c(2,i+1) - 2.*divdf1 c(3,i) = (divdf1 - c(2,i) - divdf3)/dtau 20 c(4,i) = (divdf3/dtau)/dtau ! ! estimate max.interpolation error on (-1,1). errmax = 0.

3

Dan RacoŃi Elemente de Fortran 95 do 30 i=2,n dx = (tau(i)-tau(i-1))/step do 30 j=1,istep h = float(j)*dx evaluate (i-1)st cubic piece

! !

pnatx = c(1,i-1)+h*(c(2,i-1)+h*(c(3,i-1)+h*c(4,i-1))) ! 30

errmax = amax1(errmax,abs(g(tau(i-1)+h)-pnatx)) aloger = alog(errmax) if (n .gt. 2) decay = & (aloger - algerp)/alog(float(n)/float(n-2)) algerp = aloger 40 print 640,n,errmax,decay 640 format(i3,e12.4,f11.2) stop end

1.0 Constante şi variabile In limbajul FORTRAN sunt definite următoarele tipuri de constante: a) Constante de tip întreg, numere cu semn fără punct zecimal: 12, +54, -321 b) Constante de tip real simplă precizie, numere cu semn şi punct zecimal şi eventual exponent: 1.1 , +5.3, 8.31451e0, 9.7654e-3 c) Constante de tip real dublă precizie, numere cu semn şi punct zecimal. Exponentul se notează cu litera d (dublă precizie): 3.14d0, 8.31451d0,7.77d+3,6.754d-9. d) Constante complexe simplă precizie: (parte_reală_simplă_precizie,parte_imaginară_simplă_precizie) (1.2,-3.4) e) Constante complexe dublă precizie: (parte_reală_dublă_precizie,parte_imaginară_dublă_precizie) (+3,1234d2,+0.987d0) f) Constante de tip logic, .true. pentru adevărat şi .false. pentru fals. g) Constante de tip şir de caractere: ‘sir1’,’sirul_doi_de_caractere’,”Van’t Hoff” program constante implicit none integer i real a real(4) x real(8) s double precision f logical adev,fals character(4) sir1,sir2 complex z1 complex(8) z2 i=-3 a=1.1 x=3.5e-4 s=8.31451d0 f=0.123456789e-5 adev=.true.

4

Dan RacoŃi Elemente de Fortran 95 fals=.false. sir1='abcd' sir2="ab'd" z1=(4.,-3.) z2=(1.234d-2,0.98765d2) namelist/lista/i,a,x,s,f,adev,fals,sir1,sir2,z1,z2 write(*,lista) end program constante

Variabilele în FORTRAN sunt de tipul: a) Intregi integer listă_variabile integer(4) listă_variabile integer indice1,indice2,alfa integer(4) indice3,indice4 b) Reale simplă precizie: real listă_variabile real(4) listă_variabile real j,i,s1,beta real(4) gama,lambda_1 c) Reale dublă precizie: real(8) listă_variabile double precision listă_variabile real(8) w1,w2 double precision k1,h2 d) Complexe simplă precizie: complex listă_variabile complex(4) listă_variabile complex y2,z2,a1,b2 e) Complexe dublă precizie complex(8) lista_variabile double complex lista_variabile complex(8) zz1,a_y double complex ew,er f) Logice

5

Dan RacoŃi Elemente de Fortran 95 logical listă_variabile logical predicat1,predicat2 Variabilele de tip logic iau valorile adevărat .true. , şi fals .false. h) Sir de caractere character(lungime_şir) lista_variabile character(3) sir1,sir2 character(16) nume,prenume program constante1 ! declaratii explicite de variabile scalare implicit none ! se suspenda conventia implicita integer(kind=2) integer(2) integer(kind=4) integer(4) :: integer real(kind=4) real(4) real(4) real*4

c1,c2,c3,c4 d1,d2,d3,d4 e1,e2,e3,e4 b1,b2,b3,b4 a1,a2,a3,a4

i1,i2,i3,i4 j1,j2,j3,j4 k1,k2,k3,k4 ; l1,l2,l3,l4

; ; ; ; ;

integer*2 integer*2 integer*4 integer integer*4

c d e b a

; real*4 i ; real j real k ; real*4 l

real(kind=8) :: di1,di2,di3,di4 ; real*8 di real(8) dj1,dj2,dj3,dj4 ; real(8) dj real*8 dk1,dk2,dk3,dk4 ; real*8 dk real(8) :: de1,de2,de3,de4,de double precision dl1,dl2,dl3,dl4,dl complex(kind=4) z1,z2,z3,z4 complex(4) w1,w2,w3,w4 ; complex*8 w complex :: v1,v2,v3,v4 complex(kind=8) complex(8) double complex ::

dz1,dz2,dz3,dz4 dw1,dw2,dw3,dw4 dv1,dv2,dv3,dv4

logical(kind=2) logical(4) logical

q1,q2,q3,q4 r1,r2,r3,r4 p1,p2,p3,p4

c1=1 ; c2=2 ; c3=3 ; c4=4_2 c=c1+c2+c3*c4 d1=-1 ; d2=-2 ; d3=-3 ; d4=-4_2 d=d1+d2+d3-d4 e1=1201 ; e2=1202 ; e3=1203 ; e4=1204 e=e1+e2+e3+e4 a1=-1411 ; a2=-1412 ; a3=-1413_4 ; a4=-1414_4 a=a1*a2/a3+a4 namelist/lista1/c1,c2,c3,c4,c,d1,d2,d3,d4,d write(*,lista1) ; pause 1 namelist/lista2/e1,e2,e3,e4,e,a1,a2,a3,a4,a write(*,lista2) ; pause 2 i1=1.1001 ; i2=-1.21e-01_4 ; i3=4.12 ; i4=5.6

6

Dan RacoŃi Elemente de Fortran 95 i=i1*i2+i3/i4 j1=i1**2 ; j2=i2**2 ; j3=i3**2 ; j4=i4**2 ; j=j1+j2+j3+j4 k1=i1+j1 ; k2=i2+j2 ; k3=i3+j3 ; k4=i4/j4 ; k=k1/k2*k3/k4 l1=0._4 ; l2=0.01_4 ; l3=+0.02e+03_4 ; l4=0.03e-03 l=l1+l2+l3+l4 namelist/lista3/i1,i2,i3,i4,i,j1,j2,j3,j4,j namelist/lista4/k1,k2,k3,k4,k,l1,l2,l3,l4,l write(*,lista3) ; pause 3 write(*,lista4) ; pause 4 di1=1._8 ; di2=2._8 ; di3=3._8 ; di4=4.d0 di=di1/di2+di3**di4 dj1=-0.00123_8 ; dj2=.005_8 ; dj3=+06d-3 ; dj4=0.d0 dj=abs(dj1)+abs(dj2)+abs(dj3)+abs(dj4) dk1=dble(k1) dk2=dble(k2) dk3=dble(k3) dk4=dble(k4) dk=dk1+dk2+dk3+dk4 de1=dk1*dk2 ; de2=dk2-dk3 ; de3=(dk1-dk2)/(dk3+dk4) de4=de1*de2+de3 de=de1-de2-de3+de4 namelist/lista5/di1,di2,di3,di4,di,dj1,dj2,dj3,dj4,dj namelist/lista6/dk1,dk2,dk3,dk4,dk,de1,de2,de3,de4,de write(*,lista5) ; pause 5 write(*,lista6) ; pause 6 z1=(1.,2.) ; z2=(-1.e03_4,7.01_4) z3=cmplx(i1,i2) ; z4=cmplx(i3,i4) w1=conjg(z1) ;w2=conjg(z2); w3=z1+z2-(z3-z4) ; w4=abs(w3) w=w1/w2+w2*w3 namelist/lista7/z1,z2,z3,z4,w1,w2,w3,w4,w write(*,lista7) ; pause 7 j1=real(w) ; j2=imag(w) write(*,*)'partea reala w=',j1 write(*,*)'partea imaginara w=',j2 write(*,*)' w=',w v1=z1+w1 ; v2=z2-w2 ; v3=z3/w3 ; v4=z4*w4 namelist/lista8/v1,v2,v3,v4 write(*,lista8) ; pause 8 q1=.true. ; q2=.false. q3=c1c4 r1=dk10.' elseif(predicat2)then write(*,*)'conditia 2 este adevarata,(a9.' else write(*,*)'Toate propozitiile logice sunt false' endif end program if4

2.1.4 Compararea şirurilor de caractere Şirurile de caractere se compară lexical definind operatorii: variabila_logică=lle(sir1,sir2) mai mic sau egal lexical variabila_logică=llt(sir1,sir2) mai mic lexical variabila_logică=lgt(sir1,sir2) mai mare egal lexical variabila_logică=lge(sir1,sir2) mai mare sau egal lexical logical l1,l2 l1=lle(‘abc’,’abd’) return .true.

2.1.5 InstrucŃiunea select case InstrucŃiunea select case are următoarea sintaxă: select case(expr) case(case_value_1) block_1 case(case_value_2) block_2 . . . case default block_default end select Se evaluează expresia (expr). Rezultatul, case_index se compară cu fiecare case_value pentru a găsi o concordanŃă (poate fi numai una). Când o concordanŃă există se execută blocul respectiv şi se sare la instrucŃiunea de după end select. Dacă expresia este logică atunci: case_index .eqv. case_value Dacă expresia este de tip numeric sau de tip caracter atunci: case_index = = case_value Exemplu: integer i select case(i) case(:-1) ! i= 1

Exemplu: select case (itest.eq.1) case(.true.) call subprogram1 case(.false.) call subprogram2 end select Exemplu: select case(itest) case(1) call sub1 case(2) call sub2 case default call subd end select

Interval low: :high low:high

A match case case_index>=low case_index