Laynetworks  
Web laynetworks.com Google
Home | Site Map | Tell a friends
Management Tutorials
Download
Tutorials
History
Computer Science
Networking
OS - Linux and Unix
Source Code
Script & Languages
Protocols
Glossary
IGNOU
Quiz
About Us
Contact Us
Feedback
 
Sign up for our Email Newsletter
 
Get Paid for Your Tech Turorials / Tips

 
Home > Computer Science > Numerical and Statistical Computing
 
CS 01 CS 02 CS 03 CS 04 CS 05 CS 06 CS 07 CS 08 CS 09 CS 10 CS 11 CS 12 CS 13 CS 14 CS 15 CS 16 CS 17
Page : 1 2 3 4 5 6 7
Numerical and Statistical Computing
 

PROJECT

(c) Write a program for each of the above methods using FORTRAN 90 as programming language. - (7 Marks)

Fortran 90 code for solving linear equations using Gauss Elimination without pivoting technique, gauss.f90

 

! ---------------------PROGRAM GAUSSIAN

! gauss.f90

! PROGRAM OF GAUSSIAN ELIMINATION

PARAMETER (IN=20)

REAL:: A(IN,IN), B(IN)

PRINT *,'INPUT NUMBER OF EQUATIONS N='

READ (5,*) N

PRINT *,'INPUT MATRIX COEFFICIENTS A(I,J)='

READ (5,*) ((A(I,J),J=1,N),I=1,N)

PRINT *,'INPUT RIGHT-HAND SIDE VECTOR B(I)='

READ (5,*) (B(I),I=1,N)

WRITE (6,*) ('****GAUSSIAN ELIMINATION ****')

WRITE (6,*)

WRITE (6,*) ('COEFFICIENT MATRIX you inputed:')

CALL PRINTA(A,IN,N,N,6)

WRITE(6,*)

WRITE(6,*) ('right hand side vector you inputed:')

CALL PRINTV(B,N,6)

WRITE(6,*)

! CONVERT TO UPPER TRIANGULAR FORM

DO K = 1,N-1

IF (ABS(A(K,K)).GT.1.E-6) THEN

DO I = K+1, N

X = A(I,K)/A(K,K)

DO J = K+1, N

A(I,J) = A(I,J) -A(K,J)*X

ENDDO

B(I) = B(I) - B(K)*X

ENDDO

ELSE

WRITE (6,*) 'ZERO PIVOT FOUND IN LINE:'

WRITE (6,*) K

STOP

END IF

ENDDO

WRITE(6,*) 'MODIFIED MATRIX'

CALL PRINTA(A,IN,N,N,6)

WRITE(6,*)

WRITE(6,*) 'MODIFIED RIGHT HAND SIDE VECTOR'

CALL PRINTV (B,N,6)

WRITE(6,*)

! BACK SUBSTITUTION

DO I = N,1,-1

SUM = B(I)

IF (I.LT.N) THEN

DO J= I+1,N

SUM = SUM - A(I,J)*B(J)

ENDDO

END IF

B(I) = SUM/A(I,I)

ENDDO

! PRINT THE RESULTS

write(6,*) ('SOLUTION VECTOR')

CALL PRINTV(B,N,6)

!

END PROGRAM GAUSSIAN

!------------------------------------------

 

SUBROUTINE PRINTA(A,IA,M,N,ICH)

!

! WRITE A 2D ARRAY TO OUTPUT CHANNEL 'ICH'

!

REAL A(IA,*)

DO I =1,M

WRITE(ICH,2) (A(I,J),J=1,N)

ENDDO

2 FORMAT(1X,6E12.4)

!

END SUBROUTINE PRINTA

!-----------------------------------------

 

SUBROUTINE PRINTV(VEC,N,ICH)

!

! WRITE A COLUMN VECTOR TO CHANNEL 'ICH'

!

REAL VEC(*)

WRITE(ICH,1) (VEC(I),I=1,N)

1 FORMAT(1X,6E12.4)

!

END SUBROUTINE PRINTV

!-----------------------------------------

Compiling, input data, and executing the objected code
%f90 -o run gauss.f90
%run
INPUT NUMBER OF EQUATIONS N=3
INPUT MATRIX COEFFICIENTS A(I,J)=
10.0 1.0 -5.0
-20.0 3.0 20.0
5.0 3.0 5.0
INPUT RIGHT-HAND SIDE VECTOR B(I)=1.0 2.0 6.0

The sample of final results can be displayed as

**** GAUSSIAN ELIMINATION ****

Coefficient matrix you inputed:
  .1000E+02 .1000E+01 -.5000E+01
-.2000E+02 .3000E+01 .2000E+02
.5000E+01 .3000E+01 .5000E+01
right hand side vector you inputed:
  .1000E+01 .2000E+01 .6000E+01

Modified Matrix:

  .1000E+02 .1000E+01 -.5000E+01
-.2000E+02 .5000E+01 .1000E+02
.5000E+01 .2500E+01 .2500E+01
Modified right hand side vector:
  .1000E+01 .4000E+01 .3500E+01
Solution vector:
  .1000E+01 -.2000E+01 .1400E+01

 

If you want to input data from a file rather than key strikes, you can create a data file, say, gauss.dat:

3
10.0 1.0 -5.0
-20.0 3.0 20.0
5.0 3.0 5.0
1.0 2.0 6.0

Then you can use the following commands.

%f90 -o run gauss.f90
%run<gauss.dat

Fortran 90 code for solving linear equations using Gauss-Seidel iterative method

 

PROGRAM GAUSS_SEIDEL

! Gauss-Seidel method for linear algebraic euqations

! gauss_seidel.f90

!

PARAMETER (IN=50)

REAL A(IN,IN), B(IN), X(IN), XNEW(IN), U(IN,IN)

PRINT *,'INPUT NUMBER OF EQUATIONS N='

READ (5,*) N

PRINT *,'INPUT MATRIX COEFFICIENTS A(I,J)='

READ (5,*) ((A(I,J),J=1,N),I=1,N)

PRINT *,'INPUT RIGHT-HAND SIDE VECTOR B(I)='

READ (5,*) (B(I),I=1,N)

WRITE(*,*) 'Input your gessed values:'

READ(5,*) (X(I),I=1,N)

PRINT *,'Input your tolence and max iteration:'

READ(5,*) TOL,ITS

WRITE(6,*) ('******Gauss-Seidel Method******')

WRITE(6,*)

WRITE(6,*) ('COEFFICIENT MATRIX')

WRITE(6,*)

CALL PRINTA(A,IN,N,N,6)

WRITE(6,*)

WRITE(6,*) ('RIGHT HAND VECTOR')

WRITE(6,*)

CALL PRINTV(B,N,6)

WRITE(6,*)

 

DO I=1,N

DIAG=A(I,I)

DO J=1,N

A(I,J)=A(I,J)/DIAG

ENDDO

B(I)=B(I)/DIAG

ENDDO

CALL NULL(U,IN,N,N)

DO I=1,N

DO J=I+1,N

U(I,J)=-A(I,J)

A(I,J)=0.0

ENDDO

ENDDO


WRITE(6,*) ('FIRST FEW ITERATIONS')

ITERS=0

4 ITERS=ITERS+1

CALL MVMULT(U,IN,X,N,N,XNEW)

CALL VECADD(B,XNEW,XNEW,N)

CALL SUBFOR(A,IN,XNEW,N)

CALL CHECON(XNEW,X,N,TOL,ICON)

IF (ITERS.LE.5) CALL PRINTV(X,N,6)

IF (ICON.EQ.0.AND.ITERS.LT.ITS) GO TO 4

WRITE(6,*)

WRITE(6,*) ('ITERATIONS TO CONVERGENCE')

WRITE(6,*) ITERS

WRITE(6,*)

WRITE(6,*) ('SOLUTION VECTOR')

CALL PRINTV(X,N,6)

END PROGRAM GAUSS_SEIDEL

 

SUBROUTINE CHECON(LOADS,OLDLDS,N,TOL,ICON)

REAL LOADS(*),OLDLDS(*)

ICON=1

BIG=0.

DO I=1,N

IF (ABS(LOADS(I)).GT.BIG) BIG=ABS(LOADS(I))

ENDDO

DO I=1,N

IF (ABS(LOADS(I)-OLDLDS(I))/BIG.GT.TOL) ICON=0

OLDLDS(I) = LOADS(I)

ENDDO

END SUBROUTINE CHECON

 

SUBROUTINE NULL(A,IA,M,N)

REAL A(IA,*)

DO I=1,M

DO J=1,N

A(I,J)=0.0

ENDDO

ENDDO

END SUBROUTINE NULL

 

SUBROUTINE PRINTA(A,IA,M,N,ICH)

REAL A(IA,*)

DO I=1,M


WRITE (ICH,2) (A(I,J),J=1,N)

ENDDO

2 FORMAT (1X,6E12.4)

END SUBROUTINE PRINTA

SUBROUTINE PRINTV(VEC,N,ICH)


REAL VEC(*)

WRITE(ICH,1) (VEC(I),I=1,N)

1 FORMAT (1X,6E12.4)

END SUBROUTINE PRINTV

 

SUBROUTINE MVMULT(M,IM,V,K,L,Y)

REAL M(IM,*),V(*),Y(*)

DO I=1,K

X=0.

DO J=1,L

X=X+M(I,J)*V(J)

ENDDO

Y(I)=X

ENDDO

END SUBROUTINE MVMULT

 

SUBROUTINE VECADD(A,B,C,N)

REAL A(*),B(*),C(*)

DO I=1,N

C(I)=A(I)+B(I)

ENDDO

END SUBROUTINE VECADD

 

SUBROUTINE SUBFOR(A,IA,B,N)

REAL A(IA,*), B(*)

DO I=1,N

SUM=B(I)

IF (I.GT.1) THEN

DO J=1,I-1

SUM=SUM-A(I,J)*B(J)

ENDDO

END IF

B(I)=SUM/A(I,I)

ENDDO

END SUBROUTINE SUBFOR

 

 

SUBROUTINE MSMULT(A,IA,C,M,N)

REAL A(IA,*)

DO I=1,N

DO J=1,N

A(I,J)=A(I,J)*C

ENDDO

ENDDO

END SUBROUTINE MSMULT


Compiling, input data, and executing the objected code


%f90 -o run gauss_seidel.f90

%run

INPUT NUMBER OF EQUATIONS N=

3

INPUT MATRIX COEFFICIENTS A(I,J)=

16 4 8

4 5 -4

8 -4 22

INPUT RIGHT-HAND SIDE VECTOR B(I)=

4 2 5

Input your gessed values:

1.0 1.0 1.0

Input your tolence and max iteration:

1.e-5

100

The sample of final results can be displayed as


******Gauss-Seidel Method******

 

COEFFICIENT MATRIX

.1600E+02 .4000E+01 .8000E+01

.4000E+01 .5000E+01 -.4000E+01

.8000E+01 -.4000E+01 .2200E+02


RIGHT HAND VECTOR

.4000E+01 .2000E+01 .5000E+01


FIRST FEW ITERATIONS

-.5000E+00 .1600E+01 .7000E+00

-.5000E+00 .1360E+01 .6564E+00

-.4182E+00 .1260E+01 .6084E+00

-.3691E+00 .1182E+01 .5764E+00

-.3337E+00 .1128E+01 .5537E+00


ITERATIONS TO CONVERGENCE

30


SOLUTION VECTOR

-.2500E+00 .1000E+01 .5000E+00

 

TOP
 
Page : 1 2 3 4 5 6 7
   
Donation | Useful links | Link to Laynetworks.com | Legal
Copyright © 2000-2010 Lay Networks All rights reserved.