TDMA 11.f90 - Solution of system of linear equatrions by Thomas method
From CFD-Wiki
(Difference between revisions)
Line 1: | Line 1: | ||
+ | <pre> | ||
+ | |||
+ | !Sample program for solving Smith-Hutton Test using different schemes | ||
+ | !of covective terms approximation - TDMA algorithm modul | ||
+ | !Copyright (C) 2005 Michail Kirichkov | ||
+ | |||
+ | !This program is free software; you can redistribute it and/or | ||
+ | !modify it under the terms of the GNU General Public License | ||
+ | !as published by the Free Software Foundation; either version 2 | ||
+ | !of the License, or (at your option) any later version. | ||
+ | |||
+ | !This program is distributed in the hope that it will be useful, | ||
+ | !but WITHOUT ANY WARRANTY; without even the implied warranty of | ||
+ | !MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the | ||
+ | !GNU General Public License for more details. | ||
+ | |||
+ | !You should have received a copy of the GNU General Public License | ||
+ | !along with this program; if not, write to the Free Software | ||
+ | !Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. | ||
+ | |||
+ | !*************************************************************************** | ||
+ | |||
Subroutine TDMA_1(NF) | Subroutine TDMA_1(NF) | ||
Line 6: | Line 28: | ||
!-------------------------------------------------------------------------------------------------------------------------------------------------------- | !-------------------------------------------------------------------------------------------------------------------------------------------------------- | ||
- | |||
!-------------------------------------------------------------------------------------------------------------------------------------------------------- | !-------------------------------------------------------------------------------------------------------------------------------------------------------- | ||
- | |||
!-------------------------------------------------------------------------------------------------------------------------------------------------------- | !-------------------------------------------------------------------------------------------------------------------------------------------------------- | ||
Do 101 J = 2, NYmax | Do 101 J = 2, NYmax | ||
- | + | P(1) = 0. | |
- | + | ||
- | + | ||
Q(1) = F(1,j,nf) | Q(1) = F(1,j,nf) | ||
- | |||
P(NXmaxP) = 0. | P(NXmaxP) = 0. | ||
- | |||
Q(NXmaxP) = F(NXmaxP,j,nf) | Q(NXmaxP) = F(NXmaxP,j,nf) | ||
- | |||
! Forward Elimination | ! Forward Elimination | ||
Do 10 i = 2,NXmaxP-1 | Do 10 i = 2,NXmaxP-1 | ||
- | + | temp = Ap(i,j) - Aw(i,j) * P(i-1) | |
- | + | Spp= Sp(i,j) + As(i,j) * F(i,j-1,nf) + & | |
- | + | An(i,j) * F(i,j+1,nf) | |
P(i) = Ae(i,j) / temp | P(i) = Ae(i,j) / temp | ||
- | + | Q(i) = (Spp + Aw(i,j)*Q(i-1)) / temp | |
10 continue | 10 continue | ||
Line 50: | Line 65: | ||
Do 301 J = NYmax,2,-1 | Do 301 J = NYmax,2,-1 | ||
- | + | P(1) = 0. | |
- | + | ||
Q(1) = F(1,j,nf) | Q(1) = F(1,j,nf) | ||
- | |||
P(NXmaxP) = 0. | P(NXmaxP) = 0. | ||
- | |||
Q(NXmaxP) = F(NXmaxP,j,nf) | Q(NXmaxP) = F(NXmaxP,j,nf) | ||
Line 62: | Line 74: | ||
Do 32 i = 2,NXmaxP-1 | Do 32 i = 2,NXmaxP-1 | ||
- | + | temp = Ap(i,j) - Aw(i,j) * P(i-1) | |
- | + | Spp= Sp(i,j) + As(i,j) * F(i,j-1,nf) + & | |
- | + | An(i,j) * F(i,j+1,nf) | |
P(i) = Ae(i,j) / temp | P(i) = Ae(i,j) / temp | ||
- | + | Q(i) = (Spp + Aw(i,j)*Q(i-1)) / temp | |
- | + | ||
- | + | ||
32 continue | 32 continue | ||
Line 89: | Line 99: | ||
Do 201 I = NXmax,2,-1 | Do 201 I = NXmax,2,-1 | ||
- | + | P(1) = 0. | |
Q(1) = F(i,1,nf) | Q(1) = F(i,1,nf) | ||
- | |||
- | |||
P(NYmaxP) = 0. | P(NYmaxP) = 0. | ||
- | |||
Q(NYmaxP) = F(i,NYmaxP,nf) | Q(NYmaxP) = F(i,NYmaxP,nf) | ||
Line 101: | Line 108: | ||
Do 14 j = 2,NYmaxP-1 | Do 14 j = 2,NYmaxP-1 | ||
- | + | temp = Ap(i,j) - As(i,j) * P(j-1) | |
- | + | Spp= Sp(i,j) + Aw(i,j) * F(i-1,j,nf) + & | |
- | + | Ae(i,j) * F(i+1,j,nf) | |
P(j) = An(i,j) / temp | P(j) = An(i,j) / temp | ||
- | + | Q(j) = (Spp - As(i,j)*Q(j-1)) / temp | |
- | + | ||
- | + | ||
14 continue | 14 continue | ||
Line 128: | Line 133: | ||
Do 211 I = 2,NXmax | Do 211 I = 2,NXmax | ||
- | + | P(1) = 0. | |
- | + | ||
Q(1) = F(i,1,nf) | Q(1) = F(i,1,nf) | ||
- | |||
P(NYmaxP) = 0. | P(NYmaxP) = 0. | ||
- | |||
Q(NYmaxP) = F(i,NYmaxP,nf) | Q(NYmaxP) = F(i,NYmaxP,nf) | ||
Line 140: | Line 142: | ||
Do 15 j = 2,NYmaxP-1 | Do 15 j = 2,NYmaxP-1 | ||
- | + | temp = Ap(i,j) - As(i,j) * P(j-1) | |
- | + | ||
- | + | Spp= Sp(i,j) + Aw(i,j) * F(i-1,j,nf) + & | |
- | + | Ae(i,j) * F(i+1,j,nf) | |
- | + | P(j) = An(i,j) / temp | |
- | + | Q(j) = (Spp - As(i,j)*Q(j-1)) / temp | |
15 continue | 15 continue | ||
Line 164: | Line 166: | ||
Return | Return | ||
- | + | End</pre> | |
- | End | + |
Latest revision as of 15:45, 21 September 2005
!Sample program for solving Smith-Hutton Test using different schemes !of covective terms approximation - TDMA algorithm modul !Copyright (C) 2005 Michail Kirichkov !This program is free software; you can redistribute it and/or !modify it under the terms of the GNU General Public License !as published by the Free Software Foundation; either version 2 !of the License, or (at your option) any later version. !This program is distributed in the hope that it will be useful, !but WITHOUT ANY WARRANTY; without even the implied warranty of !MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !GNU General Public License for more details. !You should have received a copy of the GNU General Public License !along with this program; if not, write to the Free Software !Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. !*************************************************************************** Subroutine TDMA_1(NF) include 'icomm_1.f90' Dimension P(nx),Q(nx) !-------------------------------------------------------------------------------------------------------------------------------------------------------- !-------------------------------------------------------------------------------------------------------------------------------------------------------- !-------------------------------------------------------------------------------------------------------------------------------------------------------- Do 101 J = 2, NYmax P(1) = 0. Q(1) = F(1,j,nf) P(NXmaxP) = 0. Q(NXmaxP) = F(NXmaxP,j,nf) ! Forward Elimination Do 10 i = 2,NXmaxP-1 temp = Ap(i,j) - Aw(i,j) * P(i-1) Spp= Sp(i,j) + As(i,j) * F(i,j-1,nf) + & An(i,j) * F(i,j+1,nf) P(i) = Ae(i,j) / temp Q(i) = (Spp + Aw(i,j)*Q(i-1)) / temp 10 continue ! Back Substitution Do 20 i = NXmaxP-1,1,-1 F(i,j,nf) = P(i)*F(i+1,j,nf) + Q(i) 20 continue 101 continue !-------------------------------------------------------------------------------------------------------------------------------------------------------- !-------------------------------------------------------------------------------------------------------------------------------------------------------- Do 301 J = NYmax,2,-1 P(1) = 0. Q(1) = F(1,j,nf) P(NXmaxP) = 0. Q(NXmaxP) = F(NXmaxP,j,nf) ! Forward Elimination Do 32 i = 2,NXmaxP-1 temp = Ap(i,j) - Aw(i,j) * P(i-1) Spp= Sp(i,j) + As(i,j) * F(i,j-1,nf) + & An(i,j) * F(i,j+1,nf) P(i) = Ae(i,j) / temp Q(i) = (Spp + Aw(i,j)*Q(i-1)) / temp 32 continue ! Back Substitution Do 30 i = NXmaxP-1,2,-1 F(i,j,nf) = P(i)*F(i+1,j,nf) + Q(i) 30 continue 301 continue !-------------------------------------------------------------------------------------------------------------------------------------------------------- !-------------------------------------------------------------------------------------------------------------------------------------------------------- Do 201 I = NXmax,2,-1 P(1) = 0. Q(1) = F(i,1,nf) P(NYmaxP) = 0. Q(NYmaxP) = F(i,NYmaxP,nf) ! Forward Elimination Do 14 j = 2,NYmaxP-1 temp = Ap(i,j) - As(i,j) * P(j-1) Spp= Sp(i,j) + Aw(i,j) * F(i-1,j,nf) + & Ae(i,j) * F(i+1,j,nf) P(j) = An(i,j) / temp Q(j) = (Spp - As(i,j)*Q(j-1)) / temp 14 continue ! Back Substitution Do 22 j = NYmaxP-1,2,-1 ! F(i,j,nf) = P(j)*F(i,j+1,nf) + Q(j) 22 continue 201 continue !-------------------------------------------------------------------------------------------------------------------------------------------------------- !-------------------------------------------------------------------------------------------------------------------------------------------------------- Do 211 I = 2,NXmax P(1) = 0. Q(1) = F(i,1,nf) P(NYmaxP) = 0. Q(NYmaxP) = F(i,NYmaxP,nf) ! Forward Elimination Do 15 j = 2,NYmaxP-1 temp = Ap(i,j) - As(i,j) * P(j-1) Spp= Sp(i,j) + Aw(i,j) * F(i-1,j,nf) + & Ae(i,j) * F(i+1,j,nf) P(j) = An(i,j) / temp Q(j) = (Spp - As(i,j)*Q(j-1)) / temp 15 continue ! Back Substitution Do 23 j = NYmaxP-1,2,-1 ! F(i,j,nf) = P(j)*F(i,j+1,nf) + Q(j) 23 continue 211 continue !-------------------------------------------------------------------------------------------------------------------------------------------------------- Return End