Gauss - Seidel 1

[milkypros]

The program i wrote below calculates the temparature distribution on a metal rod at 11 stations using Gauss-Seidel iterative method. How can i define in the program to stop the calculations when the temperatures T2 to T11 become: ?2,?3,?4,?5,?6,?7,?8,?9,?10=10,20,30,40,50,60,70,80,90


C __________________________________________________
C |__________________________________________________|
C T1 T2 T3 T4 T5 T6 T7 T8 T9 T10 T11
C
C L=1 m
C T(x)=(100.x)/L

dimension T(11)
! Initial values for temperatures
data T/0.,1.,1.,1.,1.,1.,1.,1.,1.,1.,100./
print *,'i T1 T2 T3 T4 T5 T6 T7 T8 T9 T10 T11'
print 10,0,T
print*, Give number of max iterations n'
read*, n
do 20 i=1,n
T(1) = 0
T(2) = (T(1)+T(3))/2
T(3) = (T(2)+T(4))/2
T(4) = (T(3)+T(5))/2
T(5) = (T(4)+T(6))/2
T(6) = (T(5)+T(7))/2
T(7) = (T(6)+T(8))/2
T(8) = (T(7)+T(9))/2
T(9) = (T(8)+T(10))/2
T(10) = (T(9)+T(11))/2
T(11) = 100
print 10,i,T
20 continue
10 format (f7.3,1X,f7.3,1X,f7.3,1X,f7.3,1X,f7.3,1X,f7.3,1X,f7.3,1X,
& f7.3,1X,f7.3,1X,f7.3,1X,f7.3,1X)
stop
end

 

[GerritGroot]

This looks like a typical homework task.

You'll have to define a logical variable, that says something about your convergence criteria and use it in a while loop.... ....won't say more

 

[mikrom]

In the meantime (without reading GerriGroot's answer) I enhanced your code:

[COLOR=#0000ff]C      __________________________________________________[/color]
[COLOR=#0000ff]C     |__________________________________________________|[/color]
[COLOR=#0000ff]C     T1     T2     T3     T4     T5     T6     T7     T8   T9     T10     T11[/color]
[COLOR=#0000ff]C[/color]
[COLOR=#0000ff]C     L=1 m[/color]
[COLOR=#0000ff]C     T(x)=(100.x)/L[/color]

      [COLOR=#2e8b57][b]dimension[/b][/color] T([COLOR=#ff00ff]11[/color]), T_old([COLOR=#ff00ff]11[/color]), R([COLOR=#ff00ff]11[/color])
       [COLOR=#0000ff]! Initial values for temperatures[/color]
      [COLOR=#2e8b57][b]data[/b][/color] T_old[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]0[/color].[COLOR=#804040][b]/[/b][/color]
      [COLOR=#2e8b57][b]data[/b][/color] T[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]0[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]100[/color].[COLOR=#804040][b]/[/b][/color]
      [COLOR=#804040][b]print[/b][/color] [COLOR=#804040][b]*[/b][/color],[COLOR=#ff00ff]'i  T1  T2  T3  T4  T5  T6  T7  T8  T9  T10  T11'[/color]
      [COLOR=#804040][b]print[/b][/color] [COLOR=#ff00ff]10[/color],[COLOR=#ff00ff]0[/color],T
[COLOR=#0000ff]c     accuracy[/color]
      eps [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0.00005[/color]
      [COLOR=#a020f0]call[/color] residual(T, T_old, R, [COLOR=#ff00ff]11[/color])
      i [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color]
      [COLOR=#804040][b]do[/b][/color]
[COLOR=#0000ff]c       save previous values of T to T_old[/color]
        [COLOR=#a020f0]call[/color] save_old_vec(T, T_old, [COLOR=#ff00ff]11[/color])
        T([COLOR=#ff00ff]1[/color]) [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
        T([COLOR=#ff00ff]2[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]1[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]3[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]3[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]2[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]4[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]4[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]3[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]5[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]5[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]4[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]6[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]6[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]5[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]7[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]7[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]6[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]8[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]8[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]7[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]9[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]9[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]8[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]10[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]10[/color]) [COLOR=#804040][b]=[/b][/color] (T([COLOR=#ff00ff]9[/color])[COLOR=#804040][b]+[/b][/color]T([COLOR=#ff00ff]11[/color]))[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color]
        T([COLOR=#ff00ff]11[/color]) [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]100[/color]
        [COLOR=#804040][b]print[/b][/color] [COLOR=#ff00ff]10[/color],i,T
[COLOR=#0000ff]c       compute R = T-T_old [/color]
        [COLOR=#a020f0]call[/color] residual(T, T_old, R, [COLOR=#ff00ff]11[/color])
[COLOR=#0000ff]c       exit loop if accuracy reached [/color]
        [COLOR=#804040][b]if[/b][/color] (VecNorm(R,[COLOR=#ff00ff]11[/color]) [COLOR=#804040][b]<=[/b][/color] eps) [COLOR=#804040][b]then[/b][/color]
          [COLOR=#804040][b]go to[/b][/color] [COLOR=#ff00ff]20[/color]
        [COLOR=#804040][b]end if[/b][/color]
        i [COLOR=#804040][b]=[/b][/color] i[COLOR=#804040][b]+[/b][/color][COLOR=#ff00ff]1[/color]
      [COLOR=#804040][b]end do[/b][/color]
[COLOR=#6a5acd]20[/color]    [COLOR=#804040][b]continue[/b][/color]
[COLOR=#6a5acd]10[/color]    [COLOR=#804040][b]format[/b][/color] (i5,[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],
     [highlight #ffff00][COLOR=#0000ff]&[/color][/highlight] [COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color],[COLOR=#008080]1X[/color],[COLOR=#008080]f7.3[/color])
      [COLOR=#804040][b]stop[/b][/color]
      [COLOR=#a020f0]end[/color]

[COLOR=#2e8b57][b]      real[/b][/color] [COLOR=#a020f0]function[/color] VecNorm(Vec, n)
[COLOR=#0000ff]c       computes Norm of the vector of a given dimension      [/color]
        [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
        [COLOR=#2e8b57][b]integer[/b][/color] :: n, i
[COLOR=#2e8b57][b]        real[/b][/color] :: Vec(n), s
        s [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
        [COLOR=#804040][b]do[/b][/color] i [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], n
          s [COLOR=#804040][b]=[/b][/color] s [COLOR=#804040][b]+[/b][/color] Vec(i)[COLOR=#804040][b]*[/b][/color]Vec(i)
        [COLOR=#804040][b]end do[/b][/color]
        VecNorm [COLOR=#804040][b]=[/b][/color] [COLOR=#008080]sqrt[/color](s)
      [COLOR=#a020f0]end[/color]

      [COLOR=#a020f0]subroutine[/color] save_old_vec(T, T_old, n)
[COLOR=#0000ff]c       saves T to T_old (i.e. T_old = T)[/color]
        [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
        [COLOR=#2e8b57][b]integer[/b][/color] :: n, i
[COLOR=#2e8b57][b]        real[/b][/color] :: T(n), T_old(n)
        [COLOR=#804040][b]do[/b][/color] i [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], n
          T_old(i) [COLOR=#804040][b]=[/b][/color] T(i)
        [COLOR=#804040][b]end do[/b][/color]
      [COLOR=#a020f0]end[/color]

      [COLOR=#a020f0]subroutine[/color] residual(T, T_old, R, n)
[COLOR=#0000ff]c       computes residual vector R = T-T_old[/color]
        [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
        [COLOR=#2e8b57][b]integer[/b][/color] :: n, i
[COLOR=#2e8b57][b]        real[/b][/color] :: T(n), T_old(n), R(n)
        [COLOR=#804040][b]do[/b][/color] i [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], n
          R(i) [COLOR=#804040][b]=[/b][/color] T(i) [COLOR=#804040][b]-[/b][/color] T_old(i)
        [COLOR=#804040][b]end do[/b][/color]
      [COLOR=#a020f0]end[/color]

But if it's true that this is your homework, than read something about the iterative methods, e.g.

http://en.wikipedia.org/wiki/Jacobi_method

http://en.wikipedia.org/wiki/Gauss–Seidel_method

and improve your program so, that it solves general system of linear equations given by a matrix of dimension n.

 

[mikrom]

I forgot:
You don't need to call residual before the loop

 

[milkypros]

hi
It,s not a homework. I am not even a student. I am trying to learn fortran as fast as i can by trying to solve problems i fint on the net on my spare time because i will have to use it later for my job. I am writting the basic portion of the problems but when i try to improve it to do some more, i get a lot of errors since i am new in fortran. This problem is in the book 'Engineering Analysis Interactive Methods and Programs with FORTRAN, QuickBASIC, MATLAB, and Mathematica' and its olready solved by fortran. Thanks anyway, i am trying to learn from you by studying your answers