XBeach
C:/repositories/XBeach/trunk/src/xbeachlibrary/libxbeach.F90
Go to the documentation of this file.
00001 module libxbeach_module
00002    use iso_c_binding,        only: c_int, c_char
00003    use params,               only: parameters
00004    use spaceparamsdef,       only: spacepars
00005    use timestep_module,      only: timepars
00006    use ship_module,          only: ship
00007    use vegetation_module,    only: veggie
00008 
00009 
00010    implicit none
00011    private
00012    public executestep, outputext, final, init, getversion
00013    public par, tpar, s, sglobal
00014 #ifdef USEMPI
00015    public slocal
00016 #endif
00017    save
00018 
00019    type(parameters)                     :: par
00020    type(timepars)                       :: tpar
00021    type(spacepars), pointer             :: s
00022    type(spacepars), target              :: sglobal
00023    type(ship), dimension(:), pointer    :: sh
00024 
00025    integer                              :: n,it,error
00026    real*8                               :: tbegin
00027 
00028 #ifdef USEMPI
00029    type(spacepars), target              :: slocal
00030    real*8                               :: t0,t01
00031    logical                              :: toall = .true.
00032    logical                              :: end_program
00033    integer                              :: nxbak, nybak
00034 #endif
00035 
00036 
00037    !startinit
00038 
00039    !-----------------------------------------------------------------------------!
00040    ! Initialize program                                                          !
00041    !-----------------------------------------------------------------------------!
00042 
00043 
00044 contains
00045    integer(c_int) function init()
00046       use params,               only: params_inio, all_input
00047       use spaceparams,          only: space_alloc_scalars, ranges_init
00048       use xmpi_module
00049       use initialize_module,    only: drifter_init, wave_init, sed_init, flow_init, discharge_init, hotstart_init_1, hotstart_init_2
00050       use initialize_module,    only: setbathy_init, grid_bathy
00051       use readtide_module,      only: readtide
00052       use readwind_module,      only: readwind
00053       use beachwizard_module,   only: bwinit
00054       use timestep_module,      only: timestep_init
00055       use logging_module,       only: writelog
00056       use groundwaterflow,      only: gw_init
00057       use logging_module,       only: start_logfiles, writelog_startup
00058       use means_module,         only: means_init
00059       use output_module,        only: output_init, output
00060       use ship_module,          only: ship_init
00061       use nonh_module,          only: nonh_init
00062       use vegetation_module,    only: veggie_init
00063       use paramsconst
00064       use rainfall_module
00065 #ifdef USEMPI
00066       use xmpi_module
00067       use spaceparams
00068       use logging_module,       only: writelog_mpi
00069       use params,               only: distribute_par
00070 
00071       integer, dimension(12)               :: info
00072       character(256)                       :: line
00073       integer                              :: rank,i
00074 #endif
00075 
00076       error   = 0
00077       n = 0
00078 
00079       ! setup of MPI
00080 #ifdef USEMPI
00081       s=>slocal
00082       call xmpi_initialize
00083       call xmpi_barrier(toall)
00084       t0 = MPI_Wtime()
00085 #endif
00086 
00087       ! create log files
00088       call start_logfiles(error)
00089 
00090       ! set starting time and date
00091       call cpu_time(tbegin)
00092 
00093       ! show statup message
00094       call writelog_startup()
00095 
00096       !-----------------------------------------------------------------------------!
00097       ! Initialize simulation                                                       !
00098       !-----------------------------------------------------------------------------!
00099 
00100       ! initialize time counter
00101       it      = 0
00102 
00103       ! read input from params.txt
00104       params_inio = .false.
00105       call all_input(par)
00106 
00107       ! allocate space scalars
00108       call space_alloc_scalars(sglobal)
00109       s => sglobal
00110 
00111       ! read grid and bathymetry
00112       call grid_bathy(s,par)
00113 
00114       ! distribute grid over processors
00115 #ifdef USEMPI
00116       call xmpi_determine_processor_grid(s%nx,s%ny,par%mpiboundary,par%mmpi,par%nmpi,par%cyclic,error)
00117 #if 0
00118       ! print information about the neighbours of the processes
00119       info                      = 0
00120       info(1)                   = xmpi_orank
00121       info(2)                   = xmpi_rank
00122       info(3)                   = xmpi_prow
00123       info(4)                   = xmpi_pcol
00124       info(5)                   = xmpi_left
00125       info(6)                   = xmpi_right
00126       info(7)                   = xmpi_top
00127       info(8)                   = xmpi_bot
00128       if(xmpi_isleft)  info(9)  = 1
00129       if(xmpi_isright) info(10) = 1
00130       if(xmpi_istop)   info(11) = 1
00131       if(xmpi_isbot)   info(12) = 1
00132 
00133       do i=5,8
00134          if(info(i) .eq. MPI_PROC_NULL) then
00135             info(i) = -99
00136          endif
00137       enddo
00138 
00139       call writelog("ls"," "," ranks and neigbours (-99 means: no neighbour):")
00140       call writelog("ls"," ",' ')
00141       call writelog("ls"," ","  orank rank pcol prow left right  top  bot isleft isright istop isbot")
00142 
00143       do rank = 0,xmpi_osize-1
00144          if (rank .ne. xmpi_omaster) then
00145             call xmpi_send(rank,xmpi_imaster,info)
00146             if (xmaster) then
00147                write(line,'(i5,i5,i5,i5,i5,i6,i5,i5,i7,i8,i6,i6)') info
00148             endif
00149             call writelog("ls"," ",trim(line))
00150          endif
00151       enddo
00152       call writelog("ls"," ",' ')
00153 #endif
00154       call writelog_mpi(par%mpiboundary,error)
00155 #endif
00156 
00157       ! initialize timestep
00158       call timestep_init(par, tpar)
00159 
00160       if (xmaster) then
00161 
00162          call writelog('ls','','Initializing .....')
00163       endif
00164       if (xmaster) then
00165          if(par%hotstart==1) then 
00166             call hotstart_init_1(s,par)
00167          endif
00168       endif
00169       call setbathy_init      (s,par)
00170       ! initialize physics
00171       call readtide           (s,par)
00172       call readwind           (s,par)
00173       call flow_init          (s,par)
00174       call discharge_init     (s,par)
00175       call drifter_init       (s,par)
00176       call wave_init          (s,par)
00177       call gw_init            (s,par)
00178       call rainfall_init      (s,par)
00179       ! TODO, fix ordening of arguments....
00180       call bwinit             (s)          ! works only on master process
00181 
00182       call sed_init           (s,par)
00183 
00184       call ship_init          (s,par,sh)   ! always need to call initialise in order
00185       ! to reserve memory on MPI subprocesses.
00186       ! Note: if par%ships==0 then don't allocate
00187       ! and read stuff for sh structures
00188       call veggie_init         (s,par)
00189       !
00190       if (xmaster) then
00191          if(par%hotstart==1) then 
00192             call hotstart_init_2(s,par)
00193          endif
00194       endif
00195 #ifdef USEMPI
00196       call distribute_par(par)
00197       s => slocal
00198       !
00199       ! here an hack to ensure that sglobal is populated, also on
00200       ! the not-(o)master processes, just to get valid addresses.
00201       !
00202       if (.not. xmaster .and. .not. xomaster) then
00203          !nxbak = sglobal%nx
00204          !nybak = sglobal%ny
00205          !sglobal%nx=0
00206          !sglobal%ny=0
00207          call space_alloc_arrays_dummies(sglobal,par)
00208          !sglobal%nx = nxbak
00209          !sglobal%ny = nybak
00210       endif
00211       call space_distribute_space (sglobal,s,par)
00212 #endif
00213 
00214       call ranges_init(s)
00215 
00216       ! nonh_init does not always need to be called
00217       if (par%wavemodel==WAVEMODEL_NONH) call nonh_init(s,par)
00218 
00219       ! initialize output
00220       call means_init             (sglobal,s,par     )
00221 
00222       call output_init            (sglobal,s,par,tpar)
00223 
00224 
00225       ! store first timestep
00226       ! from this point on, xomaster will hang in subroutine output
00227       ! until a broadcast .true. is received
00228       call output(sglobal,s,par,tpar)
00229       init = 0
00230    end function init
00231 
00232    integer(c_int) function outputext()
00233       use output_module,        only: output, output_error
00234       ! store first timestep
00235       call output(sglobal,s,par,tpar,.false.)
00236       if(error==0) then
00237          outputext = 0
00238       elseif(error==1) then
00239          call output_error(s,sglobal,par,tpar)
00240          outputext = 1
00241       endif
00242 
00243    end function outputext
00244    !-----------------------------------------------------------------------------!
00245    ! Start simulation                                                            !
00246    !-----------------------------------------------------------------------------!
00247 
00248    !_____________________________________________________________________________
00249 
00250    integer(c_int) function executestep(dt)
00251       use loopcounters_module,  only: execute_counter
00252       use xmpi_module,          only: xcompute
00253       use drifter_module,       only: drifter
00254       use flow_timestep_module, only: flow
00255       use boundaryconditions,   only: wave_bc, flow_bc
00256       use morphevolution,       only: bed_update, setbathy_update, transus
00257       use wave_timestep_module, only: wave
00258       use timestep_module,      only: timestep, outputtimes_update
00259       use groundwaterflow,      only: gw_bc, gwflow
00260       use beachwizard_module,   only: assim, assim_update
00261       use ship_module,          only: shipwave
00262       use vegetation_module,    only: vegatt
00263       use wetcells_module,      only: compute_wetcells
00264       use output_module,        only: log_progress
00265       use paramsconst
00266       use rainfall_module
00267       use logging_module
00268 #ifdef USEMPI
00269       use xmpi_module
00270 #endif
00271 
00272       real*8, optional :: dt
00273 
00274 #ifdef USEMPI
00275       if (execute_counter .eq. 1) then
00276          ! exclude first pass from time measurement
00277          call xmpi_barrier
00278          t01 = MPI_Wtime()
00279       endif
00280 #endif
00281       execute_counter = execute_counter + 1
00282 
00283       executestep = -1
00284 
00285       ! determine timestep
00286       if(xcompute) then
00287 
00288          ! determine this time step's wet points
00289          call compute_wetcells(s,par)
00290          !
00291          ! determine time step
00292          call timestep(s,par,tpar,it,dt=dt,ierr=error)
00293          call outputtimes_update(par, tpar)
00294 
00295          ! update log
00296          call log_progress(par)
00297 
00298          if (error==0) then
00299             !
00300             ! Boundary conditions
00301             call wave_bc        (sglobal,s,par)
00302             if (par%gwflow==1)                                      call gw_bc          (s,par)
00303             if ((par%flow==1).or.(par%wavemodel==WAVEMODEL_NONH))   call flow_bc        (s,par)
00304             !
00305             ! Compute timestep
00306             if (par%ships==1)                                       call shipwave       (s,par,sh)
00307             if (par%swave==1)                                       call wave           (s,par)
00308             if (par%vegetation==1)                                  call vegatt         (s,par)
00309             if (par%gwflow==1)                                      call gwflow         (s,par)
00310             if ((par%flow==1).or.(par%wavemodel==WAVEMODEL_NONH))   call flow           (s,par)
00311             if (par%ndrifter>0)                                     call drifter        (s,par)
00312             if (par%sedtrans==1)                                    call transus        (s,par)
00313             if (par%bchwiz>0)                                       call assim          (s,par)             ! Beach wizard
00314             !
00315             ! Bed level update
00316             if ((par%morphology==1).and.(.not. par%bchwiz==1).and.(.not. par%setbathy==1)) call bed_update(s,par)
00317             if (par%bchwiz>0)        call assim_update   (s, par)
00318             if (par%setbathy==1)     call setbathy_update(s, par)
00319          endif
00320       endif
00321 
00322       n = n + 1
00323       executestep = 0
00324    end function executestep
00325    !_____________________________________________________________________________
00326 
00327 
00328    integer(c_int) function final()
00329       use logging_module,       only: writelog_finalize
00330       use xmpi_module
00331 
00332 
00333       !-----------------------------------------------------------------------------!
00334       ! Finalize simulation                                                         !
00335       !-----------------------------------------------------------------------------!
00336 
00337 #ifdef USEMPI
00338       end_program = .true.
00339       call xmpi_send_sleep(xmpi_imaster,xmpi_omaster) ! wake up omaster
00340       call xmpi_bcast(end_program,toall)
00341       call xmpi_barrier(toall)
00342       call writelog_finalize(tbegin,n,par%t,par%nx,par%ny,t0,t01)
00343       call xmpi_finalize
00344 #else
00345       call writelog_finalize(tbegin,n,par%t,par%nx,par%ny)
00346 #endif
00347       final = 0
00348    end function final
00349 
00350    subroutine getversion(version)
00351       character(kind=c_char,len=*),intent(inout) :: version
00352 
00353       include 'version.def'
00354       include 'version.dat'
00355 
00356       version = trim(Build_Revision)
00357    end subroutine
00358 
00359 end module libxbeach_module
 All Classes Files Functions Variables Defines