/[MITgcm]/MITgcm_contrib/cg2d_bench/ini_rhs.F
ViewVC logotype

Contents of /MITgcm_contrib/cg2d_bench/ini_rhs.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (show annotations) (download)
Fri May 12 21:58:06 2006 UTC (19 years, 2 months ago) by ce107
Branch: MAIN
Initial version of CG2D benchmark code (serial and parallel) by Chris Hill

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: $

  ViewVC Help
Powered by ViewVC 1.1.22