IOD of simple linear regression Fortran 95 프로그램

IOD of simple linear regression Fortran 95 프로그램

 

자료

 

1 1 5
1 2 8
1 3 11
1 4 14
1 5 17

설명 : x0, x1, y

 

위 자료를 iod_simple_linear_regression.dat로 저장

 

lhs와 rhs를 구하는 프로그램

 

PROGRAM iod_slr_setup

! program name
! : iteration on data of simple linear regreesion
! prgrammer
! : Park Byoungho
! date
! : 2009.10.27.

IMPLICIT NONE

! data dictionary

INTEGER, PARAMETER :: no_of_rc = 2 ! number of regression coefficient
CHARACTER(LEN=256), PARAMETER :: data_filename = 'iod_simple_linear_regression.dat' ! data file name

REAL(KIND = 8), ALLOCATABLE :: dataline(:) ! each data line when reading data

REAL(KIND = 8), ALLOCATABLE :: rhs(:) ! left-hand side, right-hand side

INTEGER :: i ! for loop
INTEGER :: status ! I/O status
CHARACTER(LEN = 40) :: error_msg ! error message

ALLOCATE(dataline(no_of_rc + 1)) ! each data line when reading data
ALLOCATE(rhs(no_of_rc)) ! right-hand side

rhs = 0

! open data file
OPEN(UNIT = 10, FILE = data_filename, STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)

IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'Data file open failed -- error message : ', error_msg
STOP
END IF

! open data file for writing left-hand side
OPEN(UNIT = 20, FILE = 'lhs.dat', STATUS = 'REPLACE', ACTION = 'WRITE')

! read each line
DO
READ(UNIT = 10, FMT = *, IOSTAT = status) dataline
IF (status /= 0) EXIT ! end of file

WRITE(20,*) '1', dataline(1)*dataline(1), '2', dataline(1)*dataline(2)
rhs(1) = rhs(1) + dataline(1)*dataline(3)
WRITE(20,*) '2', dataline(2)*dataline(2), '1', dataline(2)*dataline(1)
rhs(2) = rhs(2) + dataline(2)*dataline(3)

END DO

CLOSE(11) ! close data file

! open file for right hand side
OPEN(UNIT = 30, FILE = 'rhs.dat', STATUS = 'REPLACE', ACTION = 'WRITE')

! write right-hand side
DO i = 1, no_of_rc
WRITE(30, *) rhs(i)
END DO

CLOSE(30) ! close

END PROGRAM iod_slr_setup

 

위 소스를 iod_simple_linear_regression_setup.f95로 저장

 

LHS 정렬 명령어(꼭 unix sort 사용)

 

sort -n -o sorted_lhs.dat lhs.dat

 

iteration으로 해를 구하는 프로그램

 

PROGRAM iod_slr_solve

! program name
! : iteration on data of simple linear regreesion
! : solve program
! prgrammer
! : Park Byoungho
! date
! : 2009.10.27.

IMPLICIT NONE

! data dictionary

INTEGER, PARAMETER :: no_of_rc = 2 ! number of regression coefficient

INTEGER :: no_of_lhs ! number of left hand side lines
REAL(KIND = 8), ALLOCATABLE :: lhs(:,:), rhs(:) ! left-hand side, right-hand side, x'y
REAL(KIND = 8), ALLOCATABLE :: solutions(:) ! solutions

INTEGER :: pre_eq_no ! previous equation number
INTEGER :: cur_eq_no ! current equation number

REAL(KIND = 8) :: diag_ele ! diagonal element
REAL(KIND = 8) :: temp_rhs ! temporary rihgt hand side
REAL(KIND = 8) :: pre_sol ! previous solution

INTEGER :: iteration ! iteration
REAL(KIND = 8) :: epsilon ! sum of squares of difference between old and new solutions
INTEGER :: i ! loop
INTEGER :: status ! i/o status
CHARACTER(LEN = 40) :: error_msg ! error message
INTEGER, PARAMETER :: MAX_ITER = 100 ! maximum number of iteration
REAL(KIND = 8), PARAMETER :: criteria = 1.E-12 ! criteria for stopping

ALLOCATE(rhs(no_of_rc))
ALLOCATE(solutions(no_of_rc))

solutions = 0

!open left-hand side file
OPEN(UNIT = 10, FILE = 'sorted_lhs.dat', STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)
IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'Sorted lhs file open failed -- error message : ', error_msg
STOP
END IF

! count the lines of lhs
no_of_lhs = 0
DO
READ(UNIT = 10, FMT = *, IOSTAT = status)
IF (status /= 0) EXIT ! reach end of file
no_of_lhs = no_of_lhs + 1
END DO

ALLOCATE(lhs(no_of_lhs, no_of_rc * 2))

! write(*,*) 'number of lines of left hand side = ', no_of_lhs

! store lhs to array
REWIND(10)
DO i = 1, no_of_lhs
READ(UNIT = 10, FMT = *, IOSTAT = status) lhs(i,:)
IF (status /= 0) EXIT ! reach end of file
END DO
CLOSE(11)

DO i = 1, no_of_lhs
! write(*,*) lhs(i,:)
END DO

!open right-hand side file
OPEN(UNIT = 20, FILE = 'rhs.dat', STATUS = 'OLD', ACTION = 'READ', IOSTAT = status, IOMSG = error_msg)

IF (status /= 0) THEN ! file open failed
WRITE (*,'(1X, A, A)') 'RHS file open failed -- error message : ', error_msg
STOP
END IF

! read and store right-hand side
DO i = 1, no_of_rc
READ(UNIT = 20,FMT = *) rhs(i)
END DO
close(20)

! write(*,*) 'right hand side ...'
DO i = 1, no_of_rc
! write(*,*) i, rhs(i)
END DO

DO iteration = 1, MAX_ITER

pre_eq_no = 1
diag_ele = 0.0
temp_rhs = rhs(1)
epsilon = 0.0

! start reading and processing each line of lhs
DO i = 1, no_of_lhs

cur_eq_no = lhs(i,1)

IF(pre_eq_no /= cur_eq_no) THEN
pre_sol = solutions(pre_eq_no) ! store old solution
solutions(pre_eq_no) = temp_rhs / diag_ele ! get a new solution
epsilon = epsilon + (pre_sol - solutions(pre_eq_no)) ** 2 ! calculate sum of squares of difference between old solution and new solution

diag_ele = 0.0
temp_rhs = rhs(INT(lhs(i,1)))
END IF

diag_ele = diag_ele + lhs(i, 2)
temp_rhs = temp_rhs - lhs(i, 4) * solutions(INT(lhs(i,3)))
pre_eq_no = lhs(i,1)

END DO
! end reading lhs

! calculate last solution
pre_sol = solutions(pre_eq_no) ! store old solution
solutions(pre_eq_no) = temp_rhs / diag_ele ! get a new solution
epsilon = epsilon + (pre_sol - solutions(pre_eq_no)) ** 2 ! calculate sum of squares of difference between old solution and new solution

! WRITE(*,*) iteration,'th iteration''s solutions'
! DO i = 1, no_of_rc
! WRITE(*,*) i, solutions(i)
! END DO

epsilon = epsilon / no_of_rc

! write iteration number and epsilon
write(*,*) 'iteration = ', iteration , ', epsilon = ', epsilon
IF (epsilon < criteria) THEN
EXIT
END IF

END DO
! end iteration

! open file for writing solution
OPEN(UNIT=30, FILE='sol.dat', STATUS='REPLACE', ACTION='WRITE', IOSTAT=status)
DO i = 1, no_of_rc
WRITE(30,*) i, solutions(i)
END DO
CLOSE(30)

END PROGRAM iod_slr_solve

 

위 소스를 iod_simple_linear_regression_solve.f95로 저장

 

소스 컴파일 및 실행 화면

 

자료 및 프로그램

 

 

1256637251_iod_simple_linear_regression.zip

다운로드

1256637251_iod_simple_linear_regression.zip