/[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.2 - (show annotations) (download)
Fri May 12 22:34:55 2006 UTC (17 years, 11 months ago) by ce107
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +1 -1 lines
Added CVS Id tag

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

  ViewVC Help
Powered by ViewVC 1.1.22