source: CONFIG_DEVT/IPSLCM6.5_work_ENSEMBLES/oasis3-mct/pyoasis/examples/4-orange/fortran/receiver.F90 @ 5725

Last change on this file since 5725 was 5725, checked in by aclsce, 3 years ago

Added new oasis3-MCT version to be used to handle ensembles simulations with XIOS.

File size: 1.9 KB
Line 
1program receiver
2   use mod_oasis
3   implicit none
4   integer :: i, kinfo, date
5   integer :: comp_id, part_id, var_id
6   integer, parameter :: n_points = 16
7   integer :: part_params(OASIS_Serial_Params)
8   integer :: var_nodims(2)
9   character(len=8) :: comp_name = "receiver"
10   character(len=8) :: var_name = "FRECVATM"
11   real(kind=8) :: field(n_points), error, epsilon
12
13   print '(2A)', "Component name: ", comp_name
14
15   call oasis_init_comp(comp_id, comp_name, kinfo)
16   if(kinfo<0) call oasis_abort(comp_id, comp_name, &
17      & "Error in oasis_init_comp: ", rcode=kinfo)
18   print '(A,I0)', "Receiver: Component ID: ", comp_id
19
20   part_params(OASIS_Strategy) = OASIS_Serial
21   part_params(OASIS_Length)   = n_points
22   call oasis_def_partition(part_id, part_params, kinfo)
23   if(kinfo<0) call oasis_abort(comp_id, comp_name, &
24      & "Error in oasis_def_partition: ", rcode=kinfo)
25   print '(A,I0)', "Receiver: part_id: ", part_id
26
27   var_nodims=[1, 1]
28   print '(2A)', "Receiver: var_name: ", var_name
29   call oasis_def_var(var_id, var_name, part_id, var_nodims, OASIS_IN, &
30      &              [1], OASIS_REAL, kinfo)
31   if(kinfo<0 .or. var_id<0) call oasis_abort(comp_id, comp_name, &
32      & "Error in oasis_def_var: ", rcode=kinfo)
33   print '(A,I0)', "Receiver: var_id: ", var_id
34
35   call oasis_enddef(kinfo)
36   if(kinfo<0) call oasis_abort(comp_id, comp_name, &
37      & "Error in oasis_enddef: ", rcode=kinfo)
38
39   date=0
40
41   field(:)=0
42
43   call oasis_get(var_id, date, field, kinfo)
44   if(kinfo<0) call oasis_abort(comp_id, comp_name, &
45      & "Error in oasis_get: ", rcode=kinfo)
46
47   call oasis_terminate(kinfo)
48   if(kinfo<0) call oasis_abort(comp_id, comp_name, &
49      & "Error in oasis_terminate: ", rcode=kinfo)
50
51   epsilon=1e-8
52   error=0.
53   do i = 1, n_points
54      error=error+abs(field(i)-i)
55   end do
56   if(error<epsilon) print '(A)', "Receiver: Data received successfully"
57
58end program receiver
Note: See TracBrowser for help on using the repository browser.