| 1 | 
ce107 | 
1.1 | 
CStartOfInterface | 
| 2 | 
  | 
  | 
      SUBROUTINE INI_RHS | 
| 3 | 
  | 
  | 
C     /==========================================================\ | 
| 4 | 
  | 
  | 
C     | SUBROUTINE INI_RHS                                       | | 
| 5 | 
  | 
  | 
C     | o Initialise 2d conjugate gradient solver right-hand side| | 
| 6 | 
  | 
  | 
C     |==========================================================| | 
| 7 | 
  | 
  | 
C     | Set a source term b in Ax = b (1)                        | | 
| 8 | 
  | 
  | 
C     | We solve (1) with neumann bc's whic means that b must    | | 
| 9 | 
  | 
  | 
C     | integrate out to zero over the whole domain. If b does   | | 
| 10 | 
  | 
  | 
C     | not integrate out to zero then the solution will         | | 
| 11 | 
  | 
  | 
C     | converge, but to a non-zero final residual.              | | 
| 12 | 
  | 
  | 
C     \==========================================================/ | 
| 13 | 
  | 
  | 
      IMPLICIT NONE | 
| 14 | 
  | 
  | 
 | 
| 15 | 
  | 
  | 
C     === Global variables === | 
| 16 | 
  | 
  | 
#include "SIZE.h" | 
| 17 | 
  | 
  | 
#include "EEPARAMS.h" | 
| 18 | 
  | 
  | 
#include "PARAMS.h" | 
| 19 | 
  | 
  | 
#include "CG2D.h" | 
| 20 | 
  | 
  | 
 | 
| 21 | 
  | 
  | 
C     === Routine arguments === | 
| 22 | 
  | 
  | 
CEndOFInterface | 
| 23 | 
  | 
  | 
 | 
| 24 | 
  | 
  | 
C     === Local variables === | 
| 25 | 
  | 
  | 
C     iG, jG     - Global coordinate index | 
| 26 | 
  | 
  | 
C     faceArea   - Temporary used to hold cell face areas. | 
| 27 | 
  | 
  | 
C     I,J,K      - Loop counters | 
| 28 | 
  | 
  | 
C     iMid, jMid - Global coords of mid-point of model domain | 
| 29 | 
  | 
  | 
      INTEGER  I,  J | 
| 30 | 
  | 
  | 
      INTEGER iG, jG | 
| 31 | 
  | 
  | 
      INTEGER iMid, jMid | 
| 32 | 
  | 
  | 
      INTEGER iQ, i3Q | 
| 33 | 
  | 
  | 
 | 
| 34 | 
  | 
  | 
C--   Get model global domain mid-point | 
| 35 | 
  | 
  | 
      iMid = Nx/2 | 
| 36 | 
  | 
  | 
      iQ   = Nx/4 | 
| 37 | 
  | 
  | 
      i3Q  = 3*Nx/4 | 
| 38 | 
  | 
  | 
      jMid = Ny/2 | 
| 39 | 
  | 
  | 
 | 
| 40 | 
  | 
  | 
C--   Set dummy source term for elliptic equation | 
| 41 | 
  | 
  | 
      DO J=1-OLy,sNy+OLy | 
| 42 | 
  | 
  | 
       DO I=1-OLx,sNx+OLx | 
| 43 | 
  | 
  | 
        cg2d_b(I,J) = 0. _d 0 | 
| 44 | 
  | 
  | 
        cg2d_x(I,J) = 0. _d 0 | 
| 45 | 
  | 
  | 
       ENDDO | 
| 46 | 
  | 
  | 
      ENDDO | 
| 47 | 
  | 
  | 
      DO J=1,sNy | 
| 48 | 
  | 
  | 
       DO I=1,sNx | 
| 49 | 
  | 
  | 
        cg2d_x(I,J) = 0. _d 0 | 
| 50 | 
  | 
  | 
        cg2d_b(I,J) = 0. _d 0 | 
| 51 | 
  | 
  | 
C--     Set +/-1 source function in middle of domain. | 
| 52 | 
  | 
  | 
        iG = myXGlobalLo + I - 1 | 
| 53 | 
  | 
  | 
        jG = myYGlobalLo + J - 1 | 
| 54 | 
  | 
  | 
        IF ( iG .EQ.  1   .AND. jG .EQ. 1 ) cg2d_b(I,J) =  1. _d 0 | 
| 55 | 
  | 
  | 
        IF ( iG .EQ.  Nx  .AND. jG .EQ. Ny ) cg2d_b(I,J) = -1. _d 0 | 
| 56 | 
  | 
  | 
       ENDDO | 
| 57 | 
  | 
  | 
      ENDDO | 
| 58 | 
  | 
  | 
C--   Update overlap regions synchronously | 
| 59 | 
  | 
  | 
      CALL EXCH_XY_R8(cg2d_b) | 
| 60 | 
  | 
  | 
      CALL EXCH_XY_R8(cg2d_x) | 
| 61 | 
  | 
  | 
C | 
| 62 | 
  | 
  | 
      RETURN | 
| 63 | 
  | 
  | 
      END | 
| 64 | 
  | 
  | 
 | 
| 65 | 
  | 
  | 
C $Id: $ |