source: CPL/oasis3-mct/branches/OASIS3-MCT_2.0_branch/lib/mct/examples/simple/twocmp.seqNB.F90 @ 4775

Last change on this file since 4775 was 4775, checked in by aclsce, 4 years ago
  • Imported oasis3-mct from Cerfacs svn server (not suppotred anymore).

The version has been extracted from https://oasis3mct.cerfacs.fr/svn/branches/OASIS3-MCT_2.0_branch/oasis3-mct@1818

File size: 8.3 KB
Line 
1!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2!    Math and Computer Science Division, Argonne National Laboratory   !
3!-----------------------------------------------------------------------
4! CVS twocmp.seqNB.F90,v 1.4 2004-06-24 21:07:01 eong Exp
5! CVS MCT_2_8_0
6!BOP -------------------------------------------------------------------
7!
8! !ROUTINE:  twocmp.seqNB
9!
10! !DESCRIPTION:  Provide a simple example of using MCT to connect to
11!  components executing sequentially in a single executable using
12!  the non-blocking communications to transfer data.
13
14!
15! !INTERFACE:
16!
17      program twocmpseqNB
18!
19! !USES:
20!
21!--- Use only the things needed from MCT
22      use m_MCTWorld,only: MCTWorld_init => init
23
24      use m_GlobalSegMap,only: GlobalSegMap
25      use m_GlobalSegMap,only: MCT_GSMap_init => init
26      use m_GlobalSegMap,only: MCT_GSMap_lsize => lsize
27      use m_GlobalSegMapComms,only: MCT_GSMap_recv => recv
28      use m_GlobalSegMapComms,only: MCT_GSMap_isend => isend
29      use m_GlobalSegMapComms,only: MCT_GSMap_bcast => bcast
30
31      use m_AttrVect,only    : AttrVect
32      use m_AttrVect,only    : MCT_AtrVt_init => init
33      use m_AttrVect,only    : MCT_AtrVt_zero => zero
34      use m_AttrVect,only    : MCT_AtrVt_lsize => lsize
35      use m_AttrVect,only    : MCT_AtrVt_indexRA => indexRA
36      use m_AttrVect,only    : MCT_AtrVt_importRA => importRAttr
37
38      use m_Router,only: Router
39      use m_Router,only: MCT_Router_init => init
40
41      use m_Transfer,only : MCT_ISend => isend
42      use m_Transfer,only : MCT_Recv => recv
43
44      implicit none
45
46      include 'mpif.h'
47
48      integer,parameter :: npoints = 24  ! total number of grid points
49      integer ier,nprocs,i
50      integer color,myrank,comm1,comm2
51      integer,dimension(:),pointer :: myids
52      integer,dimension(:),pointer :: req1,req2
53!-----------------------------------------------------------------------
54!  The Main program.
55! We are implementing a single-executable, seqeuntial-execution system.
56! This small main program sets up MCTWorld, calls each "init" method
57! and then calls each component in turn.
58
59      type(GlobalSegMap) :: GSMap1,GSMap2
60      type(AttrVect) :: Av1,Av2
61
62      call MPI_init(ier)
63
64      call mpi_comm_size(MPI_COMM_WORLD, nprocs,ier)
65      call mpi_comm_rank(MPI_COMM_WORLD, myrank,ier)
66
67! Duplicate MPI_COMM_WORLD into a communicator for each model
68      call mpi_comm_dup(MPI_COMM_WORLD,comm1,ier)
69      call mpi_comm_dup(MPI_COMM_WORLD,comm2,ier)
70
71      allocate(myids(2))
72      myids(1)=1
73      myids(2)=2
74
75! Initialize MCT world
76      call MCTWorld_init(2,MPI_COMM_WORLD,comm1,myids=myids)
77
78! Initialize the models, pass in the communicators
79      call model1init(comm1,req1,GSMap1,Av1)
80      call model2init(comm2,req2,GSMap2,Av2)
81
82!-----------------end of initialization phase ------
83! Run the models, pass in the communicators
84      do i=1,5
85       write(6,*) " "
86       write(6,*) "Step ",i
87       call model1(comm1,GSMap1,Av1)
88       call model2(comm2,GSMap2,Av2)
89      enddo
90
91! Models are finished.
92      call mpi_finalize(ier)
93
94      contains
95
96!-----------------------------------------------------------------------
97!-----------------------------------------------------------------------
98! !ROUTINE:
99      subroutine model1init(comm1,req1,GSmap,av1)   ! init the first model
100
101      implicit none
102
103      integer :: comm1,mysize,ier,asize,myproc
104      integer :: fieldindx,avsize,i
105      integer,dimension(1) :: start,length
106      real,pointer :: testarray(:)
107      integer,pointer :: req1(:)
108     
109      type(GlobalSegMap) :: GSmap
110      type(AttrVect) :: av1
111!---------------------------
112
113!  find local rank and size
114      call mpi_comm_size(comm1,mysize,ier)
115      call mpi_comm_rank(comm1,myproc,ier)
116      write(6,*)myproc,"model1 size",mysize
117
118!  set up a grid and decomposition
119      asize =  npoints/mysize
120
121      start(1)= (myproc*asize) +1
122      length(1)=asize
123
124!  describe decomposition with MCT GSmap type
125      call MCT_GSMap_init(GSMap,start,length,0,comm1,1)
126
127      write(6,*)myproc,"model 1 GSMap ngseg",GSMap%ngseg,start(1)
128
129      if(myproc .eq. 0) call MCT_GSMap_Isend(GSMap,2,100,req1)
130
131!  Initialize an Attribute Vector
132      call MCT_AtrVt_init(av1,rList="field1:field2",lsize=MCT_GSMap_lsize(GSMap,comm1))
133      write(6,*)myproc,"model1 got an aV"
134
135      avsize = MCT_AtrVt_lsize(av1)
136      write(6,*)myproc,"model 1 av size", avsize
137
138      end subroutine model1init
139
140!-----------------------------------------------------------------------
141!-----------------------------------------------------------------------
142      subroutine model1(comm1,GSmap,av1)   ! run the first model
143
144      implicit none
145
146      integer :: comm1,mysize,ier,asize,myproc
147      integer :: fieldindx,avsize,i
148      integer,dimension(1) :: start,length
149      real,pointer :: testarray(:)
150     
151      type(GlobalSegMap) :: GSmap,GSmap2
152      type(AttrVect) :: av1
153      type(Router),save :: Rout
154      logical,save :: firsttime=.FALSE.
155
156      call mpi_comm_rank(comm1,myproc,ier)
157
158      if(.not.firsttime) then
159!  get other GSMap
160        if(myproc .eq. 0) call MCT_GSMap_recv(GSmap2,2,110)
161        call MCT_GSMap_bcast(GSmap2,0,comm1)
162! initialize a router
163        call MCT_Router_init(GSMap,GSmap2,comm1,Rout)
164      endif
165      firsttime=.TRUE.
166
167      avsize = MCT_AtrVt_lsize(av1)
168
169!  Fill Av with some data
170!  fill first attribute the direct way
171      fieldindx = MCT_AtrVt_indexRA(av1,"field1")
172      do i=1,avsize
173        av1%rAttr(fieldindx,i) = float(i +20*myproc)
174      enddo
175
176!  fill second attribute using Av import function
177      allocate(testarray(avsize))
178      do i=1,avsize
179        testarray(i)= cos((float(i+ 20*myproc)/npoints) * 3.14)
180      enddo
181      call MCT_AtrVt_importRA(av1,"field2",testarray)
182
183!  print out Av data
184      do i=1,avsize
185        write(6,*)myproc, "model 1 data", i,av1%rAttr(1,i),av1%rAttr(2,i)
186      enddo
187     
188!  send the data
189      call MCT_ISend(av1,Rout)
190
191
192
193      end subroutine model1
194
195!-----------------------------------------------------------------------
196!-----------------------------------------------------------------------
197! !ROUTINE:
198      subroutine model2init(comm2,req2,GSmap,av1)  ! init model 2
199
200      implicit none
201
202      integer :: comm2,mysize,ier,asize,myproc
203      integer :: i
204      integer,dimension(1) :: start,length
205      type(GlobalSegMap) :: GSmap
206      type(AttrVect) :: av1
207      integer,pointer :: req2(:)
208!---------------------------
209
210!  find local rank and size
211      call mpi_comm_size(comm2,mysize,ier)
212      call mpi_comm_rank(comm2,myproc,ier)
213      write(6,*)myproc,"model2 size",mysize
214
215!  set up a grid and decomposition
216      asize =  npoints/mysize
217
218      start(1)= (myproc*asize) +1
219      length(1)=asize
220
221!  describe decomposition with MCT GSmap type
222      call MCT_GSMap_init(GSMap,start,length,0,comm2,2)
223
224      write(6,*)myproc, "model 2 GSMap ngseg",GSMap%ngseg,start(1)
225
226      if(myproc .eq. 0) call MCT_GSMap_Isend(GSMap,1,110,req2)
227
228!  Initialize an Attribute Vector
229      call MCT_AtrVt_init(av1,rList="field1:field2",lsize=MCT_GSMap_lsize(GSMap,comm2))
230      write(6,*)myproc,"model2 got an aV"
231
232      write(6,*)myproc, "model 2 av size", MCT_AtrVt_lsize(av1)
233
234      end subroutine model2init
235
236!-----------------------------------------------------------------------
237!-----------------------------------------------------------------------
238! !ROUTINE:
239      subroutine model2(comm2,GSmap,av1)
240
241      implicit none
242
243      integer :: comm2,mysize,ier,avsize,myproc
244      integer :: i
245      integer,dimension(1) :: start,length
246      type(GlobalSegMap) :: GSmap,GSmap2
247      type(AttrVect) :: av1
248      type(Router),save   :: Rout
249      logical,save :: firsttime=.FALSE.
250!---------------------------
251
252! initialize Av to be zero everywhere
253      call MCT_AtrVt_zero(av1)
254
255      call mpi_comm_rank(comm2,myproc,ier)
256      if(.not.firsttime) then
257! receive other GSMap
258        if(myproc .eq. 0) call MCT_GSMap_recv(GSmap2,1,100)
259        call MCT_GSMap_bcast(GSmap2,0,comm2)
260!  initialize a Router
261        call MCT_Router_init(GSMap,GSmap2,comm2,Rout)
262      endif
263      firsttime=.TRUE.
264
265      avsize = MCT_AtrVt_lsize(av1)
266
267!  print out Av data before Recv
268      do i=1,avsize
269        write(6,*) myproc,"model 2 data", i,av1%rAttr(1,i),av1%rAttr(2,i)
270      enddo
271
272!  Recv the data
273      call MCT_Recv(av1,Rout)
274
275!  print out Av data after Recv.
276      do i=1,avsize
277        write(6,*) myproc,"model 2 data after", i,av1%rAttr(1,i),av1%rAttr(2,i)
278      enddo
279
280
281      end subroutine model2
282
283      end
Note: See TracBrowser for help on using the repository browser.