rbf.f90 1.2 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465
  1. program rbf
  2. ! rbf.f90
  3. ! This file creates a Radial basis function matrix.
  4. !
  5. !
  6. ! We store this matrix in matrixfile.dat
  7. !
  8. integer :: m,n,d,i,j,init(4)
  9. real(kind=8) :: l,r,coef
  10. real(kind=8), allocatable, dimension(:) :: xx,yy,value,value2
  11. real(kind=8), allocatable, dimension(:,:) :: A
  12. real ::rand
  13. open(unit=20,file='inputfile.dat',status='old')
  14. read(20,*)init
  15. n=init(1)
  16. close(20)
  17. m=n-3
  18. allocate(value(2*m),value2(m),xx(m),yy(m),A(n,n))
  19. open(unit=30,file='vecrand.dat',status='old')
  20. read(30,*)value
  21. do i=1,m
  22. xx(i)=value(2*i-1)
  23. yy(i)=value(2*i)
  24. end do
  25. close(30)
  26. do i=1,m
  27. d=0
  28. do j=1,m
  29. r=dsqrt((xx(i)-xx(j))**2+(yy(j)-yy(i))**2)
  30. if (r>0) then
  31. A(i,j)=r*r*log(r)
  32. if (d==0)then
  33. d=j
  34. coef=0.5*d*rand(0)
  35. A(i,i)=coef*coef*log(coef)
  36. end if
  37. end if
  38. end do
  39. end do
  40. do j=1,m
  41. A(m+1,j)=xx(j)
  42. A(m+2,j)=yy(j)
  43. A(m+3,j)=1.d0
  44. A(j,m+1)=xx(j)
  45. A(j,m+2)=yy(j)
  46. A(j,m+3)=1.d0
  47. end do
  48. do i=m+1,n
  49. do j=m+1,n
  50. A(i,j)=0.d0
  51. end do
  52. end do
  53. deallocate(xx,yy,value,value2)
  54. open(unit=10,file='matrixfile.dat',status='unknown')
  55. do i=1,n
  56. write(10,*) (A(i,j),j=1,n)
  57. end do
  58. close(10)
  59. deallocate(A)
  60. end program rbf