| |
| /************************************************************************ |
| * * |
| * Commonwealth Scientific and Industrial Research Organisation (CSIRO) * |
| * - Division of Information Technology (DIT) * |
| * - Division of Atmospheric Research (DAR) * |
| * * |
| * Shallow water weather model - Distributed Memory Version * |
| * * |
| * Finite difference model of shallow water equations based on :- * |
| * "The dynamics of finite difference models of the shallow water * |
| * equations" by R. Sadourney, JAS, 32, 1975. * |
| * Code from:- * |
| * "An introduction to three-dimensional climate modelling" * |
| * by Washington and Parkinson * |
| * * |
| * Programmers = David Abramson (DIT) rcoda@koel.co.rmit.oz * |
| * = Paul Whiting (DIT) rcopw@koel.co.rmit.oz * |
| * = Martin Dix (DAR) mrd@koel.co.rmit.oz * |
| * Language = BSD c using Argonne NL macros * |
| * O/S = Unix System V * |
| * H/W = Encore Multimax 320 * |
| * * |
| ************************************************************************/ |
| |
| #include <math.h> |
| #include <mpi.h> |
| #include <stdio.h> |
| #include "decs.h" |
| |
| extern void worker(); |
| MPI_Datatype * setup_res(); |
| |
| main (argc, argv) |
| int argc; |
| char * argv[]; |
| { |
| float pi=4.*(float)atan((double)1.); |
| float p[n][m]; /* Pressure (or free surface height) */ |
| float u[n][m]; /* Zonal wind */ |
| float v[n][m]; /* Meridional wind */ |
| float psi[n][m]; /* Velocity streamfunction */ |
| float pold[n][m]; |
| float uold[n][m]; |
| float vold[n][m]; |
| float h[n][m]; |
| float z[n][m]; |
| float dummy1[m]; |
| float dummy2[n][m]; |
| float tpi=pi+pi; |
| float di=tpi/(float)m; |
| float dj=tpi/(float)n; |
| int i, j, chunk_size, nxt, prv; |
| |
| int master_packet[4]; |
| float p_start[m]; |
| float u_start[m]; |
| float v_start[m]; |
| float psi_start[m]; |
| float pold_start[m]; |
| float uold_start[m]; |
| float vold_start[m]; |
| int proc_cnt; |
| int tid; |
| MPI_Datatype * res_type; |
| |
| MPI_Init(&argc, &argv); |
| MPI_Comm_size(MPI_COMM_WORLD, &proc_cnt); |
| MPI_Comm_rank(MPI_COMM_WORLD, &tid); |
| |
| if ( proc_cnt < 2 ) |
| { |
| fprintf(stderr, "must have at least 2 processes, not %d\n", proc_cnt); |
| MPI_Finalize(); |
| return 1; |
| } |
| |
| if ( (n % (proc_cnt - 1)) != 0 ) |
| { |
| if ( tid == 0 ) |
| fprintf(stderr, "(number of processes - 1) must be a multiple of %d\n", n); |
| |
| MPI_Finalize(); |
| return 1; |
| } |
| |
| if (tid != 0) { |
| worker(); |
| MPI_Barrier(MPI_COMM_WORLD); |
| MPI_Finalize(); |
| } else { |
| |
| /* master process */ |
| |
| chunk_size = n / (proc_cnt - 1); |
| |
| for (i = 1: i < proc_cnt; i++) { |
| /* calculate each worker's boundary */ |
| master_packet[JSTART] = (i - 1) * chunk_size; |
| |
| if (i == proc_cnt - 1) |
| master_packet[JEND] = n - 1; |
| else |
| master_packet[JEND] = i * chunk_size - 1; |
| |
| if (i == 1) |
| prv = proc_cnt-1; |
| else |
| prv = i-1; |
| |
| master_packet[PREV] = prv; |
| |
| if (i == proc_cnt - 1) |
| nxt = 1; |
| else |
| nxt = i+1; |
| |
| master_packet[NEXT] = nxt; |
| |
| MPI_Send(&master_packet, 4, MPI_INT, i, START_SIGNAL, |
| MPI_COMM_WORLD); |
| |
| printf("jstart=%d, jend=%d, next=%d, prev=%d\n", |
| master_packet[JSTART], |
| master_packet[JEND], |
| master_packet[NEXT], |
| master_packet[PREV]); |
| } |
| |
| |
| /* |
| initialise data structures and construct packets to be sent to workers |
| */ |
| |
| initialise(p, u, v, psi, pold, uold, vold, di, dj, z); |
| diag(1, 0., p, u, v, h, z); |
| |
| for (i = 1; i < proc_cnt; i++) { |
| for (j = 0; j < n; j++) { |
| acopy_two_to_one(p, p_start, j); |
| MPI_Send(&p_start, m, MPI_FLOAT, i, P_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(u, u_start, j); |
| MPI_Send(&u_start, m, MPI_FLOAT, i, U_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(v, v_start, j); |
| MPI_Send(&v_start, m, MPI_FLOAT, i, V_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(psi, psi_start, j); |
| MPI_Send(&psi_start, m, MPI_FLOAT, i, PSI_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(pold, pold_start, j); |
| MPI_Send(&pold_start, m, MPI_FLOAT, i, POLD_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(uold, uold_start, j); |
| MPI_Send(&uold_start, m, MPI_FLOAT, i, UOLD_ROW, |
| MPI_COMM_WORLD); |
| |
| acopy_two_to_one(vold, vold_start, j); |
| MPI_Send(&vold_start, m, MPI_FLOAT, i, VOLD_ROW, |
| MPI_COMM_WORLD); |
| } |
| } |
| |
| /* |
| receive packets back from the workers |
| */ |
| res_type = setup_res(); |
| |
| if ( debug & debug_master ) |
| printf("receiving P\n"); |
| |
| update_global_ds(res_type, P_ROW, p); |
| |
| if ( debug & debug_master ) |
| printf("receiving U\n"); |
| |
| update_global_ds(res_type, U_ROW, u); |
| |
| if ( debug & debug_master ) |
| printf("receiving V\n"); |
| |
| update_global_ds(res_type, V_ROW, v); |
| |
| if ( debug & debug_master ) |
| printf("receiving H\n"); |
| |
| update_global_ds(res_type, H_ROW, h); |
| |
| if ( debug & debug_master ) |
| printf("receiving Z\n"); |
| |
| update_global_ds(res_type, Z_ROW, z); |
| |
| for (i = 1; i < proc_cnt; i++){ |
| MPI_Send(&master_packet, 4, MPI_INT, i, END_SIGNAL, |
| MPI_COMM_WORLD); |
| } |
| |
| /* wait for all workers to end */ |
| MPI_Barrier(MPI_COMM_WORLD); |
| MPI_Finalize(); |
| |
| diag(itmax, itmax*dt, p, u, v, h, z); |
| } |
| |
| return 0; |
| } |
| |
| MPI_Datatype * |
| setup_res() |
| { |
| struct res res; |
| MPI_Aint res_disp[2]; |
| static int res_done = 0; |
| static int res_len[2] = { m, 1 }; |
| static MPI_Datatype res_old[2] = { MPI_FLOAT, MPI_INT }; |
| static MPI_Datatype res_type; |
| |
| if ( res_done ) |
| return &res_type; |
| |
| res_done++; |
| MPI_Address(&res.row[0], &res_disp[0]); |
| MPI_Address(&res.indx, &res_disp[1]); |
| res_disp[1] -= res_disp[0]; |
| res_disp[0] = 0; |
| MPI_Type_struct(2, res_len, res_disp, res_old, &res_type); |
| MPI_Type_commit(&res_type); |
| |
| return &res_type; |
| } |
| |
| /* |
| this function waits for all the workers to return the packets of |
| a particular type and then updates the master's copy of the same type |
| */ |
| update_global_ds(res_type, indx, ds) |
| MPI_Datatype * res_type; |
| int indx; |
| float ds[n][m]; |
| |
| { |
| int i; |
| int row; |
| struct res res; |
| MPI_Status status; |
| |
| for (i = 0; i < n; i++) { |
| MPI_Recv(&res, 1, *res_type, MPI_ANY_SOURCE, indx, |
| MPI_COMM_WORLD, &status); |
| |
| acopy_one_to_two(res.row, ds, res.indx); |
| } |
| |
| } |