8.4 Framework for Model and Filter on Disjoint Process Sets
8.4.3 Execution Structure of the Framework
The data assimilation for disjoint process sets for model and filter exhibits a clear separation between the model and filter parts. Both are executed concurrently on their respective processes. A flow diagram for the framework which exemplifies the SEIK filter is shown in figure 8.8. The thick green lines symbolize communication.
On execution of the program, the MPI communicators are initialized by all processes in global operations. Since in this phase of the program all processes are available, the user has to take care that the subsequent model initialization is performed only by the model processes. The allocation and initialization of model fields is not required by the filter processes. After the model initialization, the filter initialization routineFilter Init is called by all processes. In this routine, the model processes store the information on the communicators COMM MODEL and COMM COUPLE while the filter processes store the information on COMM COUPLE and COMM FILTER. Subsequently, the model processes exit the filter initialization routine. The filter processes proceed in Filter Init by allocating the arrays which are required for the chosen filter. Then the state vector x and the ensemble matrix X or the mode matrix V are initialized and sub-ensembles are distributed to all filter processes. Finally the filter processes call the filter main routine whose components are shown on the right hand side of figure 8.8.
Post−processing nsteps=0
model processes filter processes
nsteps, time, state init of fields generate mesh Initialize Model Initialize Communicators
nsteps=0 Start
Time stepper include forcing
true
include BC
Stop state
true
false false
User_Analysis
Finalize_Model Do
Next_Observation Do
nsteps>0?
Do i=1,nsteps
User_Analysis
SEIK_Analysis Filter_Init
nsteps>0?
Get_State
SEIK_Resample filterpe=1?
User_Analysis Put_State
Forecast
Figure 8.8: Flow diagram for the framework when filter and model are executed by disjoint process sets. Exemplified is the program flow for the SEIK filter. Shaded in gray are the routines of the filter framework. The thick green lines denote communi-cation. The parts of the program which are horizontally centered are executed by all processes. After the initialization, the program splits into the model part displayed on the left hand side and the filter part on the right hand side. Both parts are connected by communication operations.
Having left the routineFilter Init, the model processes proceed to the forecast loop shown on the left hand side of figure 8.8. InGet State (see algorithm 8.8), the processes wait to receive the value ofnsteps which is sent by the filter. Ifnsteps= 0, no forecast has to be performed. Thus, no further operations are necessary in Get State and the forecast loop is exited. If nsteps >0, the processes also receive the variable time and the state vector x to be evolved. Subsequently, the routine Distribute State is called which initializes the model fields on the basis of the state vectorx. Then the evolution of the state is performed by the model time stepper. After the evolution, the routine Put State is called. This routine is shown as algorithm 8.9. HereCollect State is called to initialize the forecasted state vector from the model fields on the model process with rank 0. Subsequently, this process sends the state vector x to the filter. This com-pletes the forecast loop and the processes return to the begin of the unconditioned loop.
The structure of the filter main routine on the right hand side of figure 8.8 is essen-tially the same as that of the serial algorithm which we have discussed as algorithm 3.1.
An addition to this algorithm is the subroutine Finalize Model. It is required in the parallel program to send nsteps with a value of zero to the model tasks. As discussed above, this signalizes to the model tasks to exit the forecast loop.
The subroutineForecast controls the loop over all ensemble members to be evolved.
It is shown as algorithm 8.10. In the configuration with disjoint processes for filter and model, an algorithm is used which sends a only single ensemble state vector to the available model tasks. The filter part of the algorithm uses non-blocking MPI operations. These only post the communication operation and immediately return from the function even if the communication operation is not yet completed. In contrast to this, the routines Get State and Put State apply blocking MPI operations to ensure that the data has been received or send completely. Sending and receiving single state vectors permits a flexible handling of the forecast phase. If a forecasted state vector is received back from some model task, a new ensemble state vector can be send immediately to this task if there are any ensemble states left. For sufficiently large ensembles, this ensures a good load balancing since faster model tasks can evolve more ensemble states than slower model tasks. This algorithm is more flexible than the configuration used for joint process sets for filter and model. There the sizes of sub-ensembles are set during the initialization phase of the framework. In addition, the memory requirements are smaller here. In the case of mode-decomposition, a single state vector is allocated on the model processes with rank 0 in COMM MODEL. No filter-related memory allocations are required on the remaining model processes. For domain-decomposition a single sub-state is allocated on each model process. For the configuration using joint process sets for filter and model, it is required to allocate sub-ensembles of state vectors.
Subroutine Forecast(step,nsteps,time) int step {Current time step, input}
int nsteps {Number of time steps to be computed, output} real time {Current model time, output}
int n {Model state dimension} int Np {Size of local state ensemble} real Xp(n, Np) {Local state ensemble}
int npes {Number of processes in COMM COUPLE} int status(npes−1) {Status array; idle: 0, working: 1}
int send ens {Counter for ensemble member to become evolved} int get ens {Number of received state vectors}
1: status(1 :npes−1)←0 {Set status to idle for all tasks}
2: send ens ←1 {Send first ensemble member}
3: get ens←0 {No state received yet}
4: loop
5: for task = 1, npes−1do
6: if status(task) == 1 then
7: Test whether receiving from task has been completed {With operation MPI Test}
8: if receiving oftask completedthen
9: get ens←get ens+ 1 {Increase counter of received states}
10: status(task)←0 {Set task to idle}
11: end if
12: end if
13: if status(task) == 0 then
14: send nsteps to task {With operation MPI ISend}
15: send time totask {With operation MPI ISend}
16: send Xp(:, send ens) to task {With operation MPI ISend}
17: post receiving of Xp(:, send ens) from task {Operation MPI IRecv}
18: send ens←send ens+ 1 {Increase index of member to send}
19: status(task)←1 {Set task to working}
20: end if
21: end for
22: if get ens==Np then
23: Exit loop
24: end if
25: end loop
Algorithm 8.10: Structure of the routine of the filter framework which controls the ensemble forecast in the case of SEIK and EnKF. (For SEEK, the state estimate itself is also evolved. Hence, the forecast routine for SEEK contains an extension for evolving the state estimate.) The used MPI operations are non-blocking. Thus, the algorithm directly proceeds after posting a MPI ISend or MPI IRecv operation.
} } }
0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000 0000
1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111 1111
model fields state vector
h h
}
v
u } v
}
u
Figure 8.9: Transition between the abstract state vector (left hand side) and the model