module mysubs contains subroutine init(n, x, fx, x0) use kinds use fileutil implicit none real(kind=REAL8), dimension(:), pointer :: x, fx integer, intent(out) :: n real(kind=REAL8), intent(out) :: x0 integer, parameter :: fileid = 1 character(len=256) :: fileName integer :: status, i write (*,*) 'Berechnet ein Lagrange-Polynom n-ten Grades' write (*,*) write (*,*) 'In welcher Datei liegen die Stuetzstellen?& & (In lagrange.dat liegt ein Beispiel.)' read (*,*) fileName open(unit=fileId, file=fileName, action='read', iostat=status) if (status /= 0) then write (*,*) 'Kann Datei ', fileName, 'nicht öffnen.' stop end if n = anzahlReals2(fileId) - 1 allocate(x(0:n)) allocate(fx(0:n)) read(fileid, *, iostat=status) (x(i), fx(i), i=0, n) if (status /= 0) then write (*,*) 'Konnte Array nicht lesen.' stop end if close(unit=fileid) write (*,*) 'Insgesamt ', n + 1, ' Stuetzstellen gelesen.' write (*,*) write (*,*) 'An welcher Stelle soll das Polynom berechnet werden? ' read (*,*) x0 return end subroutine init subroutine result(x, y) use kinds implicit none real(kind=REAL8), intent(in) :: x, y write (*,'(A, F10.2, A, F10.2, A)') 'Das Lagrange-Polynom an der Stelle ',& & x, ' hat den Wert ', y, '.' return end subroutine result end module mysubs program Lagrange ! Berechnet ein Lagrange-Polynom n-ten Grades anhand von n + 1 ! vorgegebenen Stuetzstellen. use kinds use mysubs implicit none real(kind=REAL8), dimension(:), pointer :: x, fx real(kind=REAL8) :: x0 ! an dieser Stelle wird das Polynom berechnet real(kind=REAL8) :: y ! enthaelt das Ergebnis real(kind=REAL8) :: produkt ! Zwischenergebnis: der i-te Produktterm integer :: n, i, j call init(n, x, fx, x0) y = 0 ! Ergebnis initialisieren do i = 0, n produkt = fx(i) do j = 0, n if (i /= j) then produkt = produkt * (x0 - x(j)) / (x(i) - x(j)) end if end do y = y + produkt end do call result(x0, y) stop end program Lagrange