examples/smpi/NAS/SP/lhsx.f

   1
   2 c---------------------------------------------------------------------
   3 c---------------------------------------------------------------------
   4
   5        subroutine lhsx(c)
   6
   7 c---------------------------------------------------------------------
   8 c---------------------------------------------------------------------
   9
  10 c---------------------------------------------------------------------
  11 c This function computes the left hand side for the three x-factors
  12 c---------------------------------------------------------------------
  13
  14        include 'header.h'
  15
  16        double precision ru1
  17        integer          i, j, k, c
  18
  19
  20 c---------------------------------------------------------------------
  21 c      treat only cell c
  22 c---------------------------------------------------------------------
  23
  24 c---------------------------------------------------------------------
  25 c      first fill the lhs for the u-eigenvalue
  26 c---------------------------------------------------------------------
  27        do  k = start(3,c), cell_size(3,c)-end(3,c)-1
  28           do  j = start(2,c), cell_size(2,c)-end(2,c)-1
  29              do  i = start(1,c)-1, cell_size(1,c)-end(1,c)
  30                 ru1 = c3c4*rho_i(i,j,k,c)
  31                 cv(i) = us(i,j,k,c)
  32                 rhon(i) = dmax1(dx2+con43*ru1,
  33      >                          dx5+c1c5*ru1,
  34      >                          dxmax+ru1,
  35      >                          dx1)
  36              end do
  37
  38              do  i = start(1,c), cell_size(1,c)-end(1,c)-1
  39                 lhs(i,j,k,1,c) =   0.0d0
  40                 lhs(i,j,k,2,c) = - dttx2 * cv(i-1) - dttx1 * rhon(i-1)
  41                 lhs(i,j,k,3,c) =   1.0d0 + c2dttx1 * rhon(i)
  42                 lhs(i,j,k,4,c) =   dttx2 * cv(i+1) - dttx1 * rhon(i+1)
  43                 lhs(i,j,k,5,c) =   0.0d0
  44              end do
  45           end do
  46        end do
  47
  48 c---------------------------------------------------------------------
  49 c      add fourth order dissipation
  50 c---------------------------------------------------------------------
  51        if (start(1,c) .gt. 0) then
  52           i = 1
  53           do   k = start(3,c), cell_size(3,c)-end(3,c)-1
  54              do   j = start(2,c), cell_size(2,c)-end(2,c)-1
  55                 lhs(i,j,k,3,c) = lhs(i,j,k,3,c) + comz5
  56                 lhs(i,j,k,4,c) = lhs(i,j,k,4,c) - comz4
  57                 lhs(i,j,k,5,c) = lhs(i,j,k,5,c) + comz1
  58
  59                 lhs(i+1,j,k,2,c) = lhs(i+1,j,k,2,c) - comz4
  60                 lhs(i+1,j,k,3,c) = lhs(i+1,j,k,3,c) + comz6
  61                 lhs(i+1,j,k,4,c) = lhs(i+1,j,k,4,c) - comz4
  62                 lhs(i+1,j,k,5,c) = lhs(i+1,j,k,5,c) + comz1
  63              end do
  64           end do
  65        endif
  66
  67        do   k = start(3,c), cell_size(3,c)-end(3,c)-1
  68           do   j = start(2,c), cell_size(2,c)-end(2,c)-1
  69              do   i=3*start(1,c), cell_size(1,c)-3*end(1,c)-1
  70                 lhs(i,j,k,1,c) = lhs(i,j,k,1,c) + comz1
  71                 lhs(i,j,k,2,c) = lhs(i,j,k,2,c) - comz4
  72                 lhs(i,j,k,3,c) = lhs(i,j,k,3,c) + comz6
  73                 lhs(i,j,k,4,c) = lhs(i,j,k,4,c) - comz4
  74                 lhs(i,j,k,5,c) = lhs(i,j,k,5,c) + comz1
  75              end do
  76           end do
  77        end do
  78
  79        if (end(1,c) .gt. 0) then
  80           i = cell_size(1,c)-3
  81           do   k = start(3,c), cell_size(3,c)-end(3,c)-1
  82              do   j = start(2,c), cell_size(2,c)-end(2,c)-1
  83                 lhs(i,j,k,1,c) = lhs(i,j,k,1,c) + comz1
  84                 lhs(i,j,k,2,c) = lhs(i,j,k,2,c) - comz4
  85                 lhs(i,j,k,3,c) = lhs(i,j,k,3,c) + comz6
  86                 lhs(i,j,k,4,c) = lhs(i,j,k,4,c) - comz4
  87
  88                 lhs(i+1,j,k,1,c) = lhs(i+1,j,k,1,c) + comz1
  89                 lhs(i+1,j,k,2,c) = lhs(i+1,j,k,2,c) - comz4
  90                 lhs(i+1,j,k,3,c) = lhs(i+1,j,k,3,c) + comz5
  91              end do
  92           end do
  93        endif
  94
  95 c---------------------------------------------------------------------
  96 c      subsequently, fill the other factors (u+c), (u-c) by a4ing to
  97 c      the first
  98 c---------------------------------------------------------------------
  99        do   k = start(3,c), cell_size(3,c)-end(3,c)-1
 100           do   j = start(2,c), cell_size(2,c)-end(2,c)-1
 101              do   i = start(1,c), cell_size(1,c)-end(1,c)-1
 102                 lhs(i,j,k,1+5,c)  = lhs(i,j,k,1,c)
 103                 lhs(i,j,k,2+5,c)  = lhs(i,j,k,2,c) -
 104      >                            dttx2 * speed(i-1,j,k,c)
 105                 lhs(i,j,k,3+5,c)  = lhs(i,j,k,3,c)
 106                 lhs(i,j,k,4+5,c)  = lhs(i,j,k,4,c) +
 107      >                            dttx2 * speed(i+1,j,k,c)
 108                 lhs(i,j,k,5+5,c) = lhs(i,j,k,5,c)
 109                 lhs(i,j,k,1+10,c) = lhs(i,j,k,1,c)
 110                 lhs(i,j,k,2+10,c) = lhs(i,j,k,2,c) +
 111      >                            dttx2 * speed(i-1,j,k,c)
 112                 lhs(i,j,k,3+10,c) = lhs(i,j,k,3,c)
 113                 lhs(i,j,k,4+10,c) = lhs(i,j,k,4,c) -
 114      >                            dttx2 * speed(i+1,j,k,c)
 115                 lhs(i,j,k,5+10,c) = lhs(i,j,k,5,c)
 116              end do
 117           end do
 118        end do
 119
 120        return
 121        end
 122
 123
 124