×
INTELLIGENT WORK FORUMS
FOR COMPUTER PROFESSIONALS

Contact US

Log In

Come Join Us!

Are you a
Computer / IT professional?
Join Tek-Tips Forums!
  • Talk With Other Members
  • Be Notified Of Responses
    To Your Posts
  • Keyword Search
  • One-Click Access To Your
    Favorite Forums
  • Automated Signatures
    On Your Posts
  • Best Of All, It's Free!

*Tek-Tips's functionality depends on members receiving e-mail. By joining you are opting in to receive e-mail.

Posting Guidelines

Promoting, selling, recruiting, coursework and thesis posting is forbidden.

Students Click Here

Gauss-Siedel Matrix to solve Elliptic Equation

Gauss-Siedel Matrix to solve Elliptic Equation

Gauss-Siedel Matrix to solve Elliptic Equation

(OP)
Qs: Write a FORTRAN program to approximately solve elliptic equation : - u_xx - u_yy=1 with the five-point finite difference formula in the right-angled triangle, sides 1,2,3^(1/2) using the Gauss-Seidel matrix solver withe the boundary conditions u=0 on all sides. Take a mesh spacing to be h=0.01. (You'll need linear interpolation on the long side.)
=>

CODE -->

PROGRAM Gauss_Seidel
IMPLICIT NONE

! Declare Variables
Real, allocatable :: x(:),y(:),u(:,:), v(:,:),u_old(:,:)
Real:: h,tolerence,error
Integer:: i,j,JI,NI

h=0.01
JI=100
NI=173
error = 1.d0
tolerence = 10E-4
 ! Total number of space stepsize
allocate (x(0:JI),y(0:NI),u(0:JI+1,0:NI+1),v(0:JI+1,0:NI+1),u_old(0:JI+1,0:NI+1))
open(10,file='Gauss_Seidel.m')  !Opening files in Matlab

!Initial Conditions
x(0)= 0
x(JI)= 1.0
y(0)= 0
y(NI)= SQRT(3.0)

do i=0,JI
do j=0,NI
x(i)= i*h   ! x-axix, x starts from 0 to 1 
y(j)= j*h   ! y-axis  y starts from 0 to SQRT(3.0)
u(i,j)= 0         ! Entire Boundary is zero
end do
end do

while (error .GT. tolerence) do  ! To stop
  do i=1, JI-1
    do j=1,NI-1
u_old(i,j)= u(i,j)  ! To store the old values

!Using 5-point scheme Formulae and rearranging the equation
v(i,j)= 0.25*(u(i+1,j)+u(i,j+1)+u(i-1,j)+u(i,j-1)+h**2) 
      end do
     end do

    do i=1, JI-1
     do j=1, NI-1
u(i,j)= v(i,j)  ! Giving the new values
        end do
      end do

error =0.d0        ! Now, error reading the value of zero
     do i=1,JI-1
       do j=1, NI-1
error = error + abs(u(i,j)- u_old(i,j))  ! To Stop
        end do
         end do
        end do

!Print out the Approximate solution in matlab to get output and plot 
write(10,*)  'x=['                           
     do i=0, JI
       write(10,*) x(i)
 end do
write(10,*) ']'

write(10,*)  'y=['                           
       do j=0,NI
write(10,*) y(j)
 end do
write(10,*) ']'


write(10,*)  'u=['                           
     do i=0, JI
       do j=0,NI
write(10,*) u(i,j)
 end do
 end do
write(10,*) ']'

write(10,*) " contour(x, y, reshape(u, length(x), length(y)))"                      !Ploting diagram x,y,u
write(10,*) "xlabel('x'),ylabel('y'),legend('Approximate Gauss Seidel')"
close(10) 

END PROGRAM Gauss_Seidel 

Note: There are no errors in the complier and no run-time errors AFTER running the program .
but my problems is h=0.01 WHICH is small and JI = 100 AND NI = 173 WHICH ARE TOO LARGE, THE FORTRAN IS REACTING THE PROGRAMS TOO SLOWY AFTER I RUN THE PROGRAM.
THE REASON WHAY I CHOOSE JI =100 AND NI= 173 IS BEACUSE X IS BETWEEN 0 TO 1 AND Y IS BETWEEN 0 TO SQUARE ROOT OF 3.
I AM REALLY CONFUSE ABOUT the five-point finite difference formula in the right-angled triangle SIDES. WHAT THESE SIDES ACTUALLY DO AND WHERE DOES IT APPLIES?
IS IT right-angled triangle SIDES GOT SOMETHING DO HERE: v(i,j)= 0.25*(u(i+1,j)+u(i,j+1)+u(i-1,j)+u(i,j-1)+h**2)
DOES MY CODE ANSWER THE QUESTION OF THE PROBLEM LIKE USING GAUSS SEIDEL MATRIX SOLVER?





RE: Gauss-Siedel Matrix to solve Elliptic Equation

(OP)

CODE

PROGRAM Gauss_Seidel
IMPLICIT NONE

! Declare Variables
Real, allocatable :: x(:),y(:),u(:,:),u_old(:,:)
Real:: h,tolerence,error
Integer:: i,j,JI,NI

h=0.01
JI=100
NI=173
error = 1.d0
tolerence = 10E-4
 ! Total number of space stepsize
allocate (x(0:JI),y(0:NI),u(0:JI+1,0:NI+1),u_old(0:JI+1,0:NI+1))
open(10,file='Gauss_Seidel.m')  !Opening files in Matlab

!Initial Conditions
x(0)= 0
x(JI)= 1.0
y(0)= 0
y(NI)= SQRT(3.0)

do i=0,JI
do j=0,NI
x(i)= i*h   ! x-axix, x starts from 0 to 1 
y(j)= j*h   ! y-axis  y starts from 0 to SQRT(3.0)
u(i,j)= 0         ! Entire Boundary is zero
end do
end do

while (error .GT. tolerence) do  ! To stop
  do i=1, JI-1
    do j=1,NI-1
u_old(i,j)= u(i,j)  ! To store the old values

!Using 5-point scheme Formulae and rearranging the equation
u(i,j)= 0.25*(u(i+1,j)+u(i,j+1)+u(i-1,j)+u(i,j-1)+h**2) 
     end do
     end do

error =0.d0        ! Now, error reading the value of zero
     do i=1,JI-1
       do j=1, NI-1
error = error + abs(u(i,j)- u_old(i,j))  ! To Stop
        end do
         end do
        end do

 do i=1, JI-1
    do j=1,NI-1

 u(i+1,j)= - (alpha * u(i,j))/ beta    ! lies outside the long side
 u(i,j+1)= - (a * u(i,j))/ b           ! lies outside the long side
 
     end do
     end do



!Print out the Approximate solution in matlab to get output and plot 
write(10,*)  'x=['                           
     do i=0, JI
       write(10,*) x(i)
 end do
write(10,*) ']'

write(10,*)  'y=['                           
       do j=0,NI
write(10,*) y(j)
 end do
write(10,*) ']'

write(10,*)  'u=['                           
     do i=0, JI
       do j=0,NI
write(10,*) u(i,j)
 end do
 end do
write(10,*) ']'

write(10,*) "[X,Y] = meshgrid(x,y)"              !Ploting diagram x,y,u
write(10,*) "Z=reshape(u,length(y), length(x))"
write(10,*) "contour(X,Y,Z)"
write(10,*) "xlabel('x'),ylabel('y'),legend('Approximate Gauss Seidel')"
close(10) 

END PROGRAM Gauss_Seidel 

THIS IS WHAT I HAVE GOT SO FAR, STILL LITTLE BIT INCOMPLETE.
THE LENGTH OF THE LONG SIDE OF THE TRIANGLE IS 2. SINCE THE POINTS ON THE LONG SIDE DO NOT COINCIDE WITH THE GRID POINTS, THE POINTS u(i,j+1) and u(i+1,j) lies outside the long side of the triangle, the distance between the points u(i,j) to u(i,j+1) and u(i,j) to u(i+1,j) are zero which is known using the linear interpolation. Alpha is the distance from zero to point u(i + 1, j) and beta is the distance from point u(i, j) to zero.
Similarly, a is the distance from zero to point u(i, j+1) and b is the distance from point u(i, j) to zero.
The points u(i -1, j) and u(i, j - 1) lies inside the long side of the triangle, the distance from these points to u(i, j) remains unchanged and does not affect in the five-point scheme.
Can anyone help me after this. What can do with alpha, beta,a and b? Do I need to know values for alpha,beta, a and b?
Please help me


Red Flag This Post

Please let us know here why this post is inappropriate. Reasons such as off-topic, duplicates, flames, illegal, vulgar, or students posting their homework.

Red Flag Submitted

Thank you for helping keep Tek-Tips Forums free from inappropriate posts.
The Tek-Tips staff will check this out and take appropriate action.

Reply To This Thread

Posting in the Tek-Tips forums is a member-only feature.

Click Here to join Tek-Tips and talk with other members! Already a Member? Login

Close Box

Join Tek-Tips® Today!

Join your peers on the Internet's largest technical computer professional community.
It's easy to join and it's free.

Here's Why Members Love Tek-Tips Forums:

Register now while it's still free!

Already a member? Close this window and log in.

Join Us             Close