Commit 33e8783d authored by Zolnierczuk, Piotr's avatar Zolnierczuk, Piotr
Browse files

added polynom test 

silenced uninitialized variables in polynom and minpack
parent 01d35734

sources/minpack.f90

+1 −0
Original line number Diff line number Diff line
@@ -270,6 +270,7 @@ 
      iflag = 0 

      nfev = 0 

      temp = zero ! (paz) initialize temp

      xnorm = zero ! (paz) initialize xnorm

!                                                                       

!     check the input parameters for errors.                            

!

sources/polynom.f90

+2 −1
Original line number Diff line number Diff line
@@ -41,7 +41,8 @@ contains 
    !

    external :: DGECO, DGESL





    am=0

    b =0

    ifail = -1

tests/Makefile

+2 −1
Original line number Diff line number Diff line
@@ -50,7 +50,8 @@ TESTSRC1=\ 
       teststrings.f90 \

       testreaddata.f90 \

	   testreduce.f90 \

	   testdrutils.f90

	   testdrutils.f90 \

	   testpolynom.f90



# user driven test programs

TESTSRC2=\

tests/testpolynom.f90

0 → 100644
+69 −0
Original line number Diff line number Diff line 
program test_polynom

  use, intrinsic :: iso_fortran_env, only : output_unit



  use polynom

  use logger



  implicit none



  double precision , parameter :: EPS=1e-8

  integer, parameter :: NDATA   =10

  integer, parameter :: MAXPOLY = 5



  double precision   :: x(NDATA), y(NDATA)

  double precision   :: a(0:MAXPOLY)

  double precision   :: aexp(0:MAXPOLY)

  double precision   :: r



  integer :: i, n



  call loginit(LOG_TRACE, prefix='test_')

  write(*,'(a)') '===> test_polynom'



  ! initialization

  x = [ (1e-2*i,i=1,NDATA) ]

  r = 0



  n = 3

  y = 4*x**3 + 3*x**2 + 2*x + 1

  aexp = [1,2,3,4,0,0]

  call test_polyfit(x, y, n, aexp(0:n))



  n = 2

  y = x**2 + 2*x + 1

  aexp = [1,2,1,0,0,0]

  call test_polyfit(x, y, n, aexp(0:n))



contains



subroutine test_polyfit(x, y, n, aexp)

    integer, intent(inout) :: n

    double precision, intent(in) :: x(:), y(:)

    double precision, intent(in) :: aexp(0:n)

      

    double precision :: yf, delta

    integer :: i

    character(len=64)  :: cmsg

        

    call polyfit(x,y,NDATA,a, n, r)

    write(output_unit,'("polyfit n=",i0)') n

    do i=n,0,-1

    delta = abs((a(i)-aexp(i))/aexp(i))

    write(cmsg,'("a[",i2,"]",g14.8," ",g14.8, " d="g14.8)') i,a(i),aexp(i),delta

    write(output_unit,*) trim(cmsg)

    call assert(delta<EPS,'polyfit '//trim(cmsg))

    end do

    write(output_unit,'("r=",g14.8)') r

    call assert(abs(r)<EPS,'polyfit-r')

    

    do i=1, NDATA

    yf = polyeval(x(i),a,n)

    delta = abs((y(i)-yf)/yf)

    write(cmsg,'("y[",i2,"]",g14.8," ",g14.8, " d="g14.8)') i, y(i),yf, delta

    write(output_unit,*) trim(cmsg)

    call assert(abs(y(i)-yf)<EPS,'polyeval '//trim(cmsg))

    end do

end subroutine test_polyfit





end program test_polynom