module libxbeach_module use iso_c_binding use params use spaceparams use xmpi_module use initialize use boundaryconditions use drifter_module use flow_timestep_module use morphevolution use readtide_module use readwind_module use wave_stationary_module use wave_timestep_module use timestep_module use readkey_module use groundwaterflow use logging_module use means_module use output_module use mnemiso_module IMPLICIT NONE ! Global variables type(parameters), save :: par real*8, save :: tbegin,tend integer, save :: it logical, save :: newstatbc type(spacepars), save, pointer :: s type(spacepars), save, target :: sglobal type(spacepars), save, target :: slocal type(timepars), save :: tpar contains integer(c_int) function init() bind(C, name="init") !DEC$ ATTRIBUTES DLLEXPORT::Init character(len=80) :: dummystring character(len=8) :: date character(len=10) :: time character(len=5) :: zone integer :: it,error #ifdef USEMPI real*8 :: t0,t01,t1 #endif character(len=155) :: cwd ! for printing the working dir ! autotools #ifdef HAVE_CONFIG_H #include "config.h" #endif ! subversion information include 'version.def' ! ---------------------------- ! Initialize program ! ---------------------------- ! Error code init = 0 include 'version.dat' ! Setup of MPI #ifdef USEMPI s=>slocal call xmpi_initialize t0 = MPI_Wtime() #endif ! Start up log files call start_logfiles(error) if (error==1) then write(*,*) 'Error: not able to open log file. Please contact XBeach team. Stopping simulation' stop endif ! call cpu_time(tbegin) call DATE_AND_TIME(DATE=date, TIME=time, ZONE=zone) ! only run this on linux (gcc function) #ifdef HAVE_CONFIG_H call getcwd(cwd) #endif if (xmaster) then call writelog('ls','','**********************************************************') call writelog('ls','',' Welcome to XBeach ') call writelog('ls','',' ') call writelog('ls','',' revision ',trim(Build_Revision) ) call writelog('ls','',' date ',trim(Build_Date) ) call writelog('ls','',' URL: ',trim(Build_URL) ) call writelog('ls','','**********************************************************') call writelog('ls','',' ') call writelog('ls','','Simulation started: YYYYMMDD hh:mm:ss time zone (UTC)') call writelog('ls','',' '//date //' '//time(1:2)//':'//time(3:4)//':'//time(5:6)//' '//zone) call writelog('ls','',' ') #ifdef HAVE_CONFIG_H call writelog('ls','',' running in: ',cwd ) #endif call writelog('ls','','General Input Module') #ifdef USEMPI if(xmaster) then call writelog('ls','','MPI version, running on ',xmpi_size,'processes') endif #endif endif ! ---------------------------- ! Initialize simulation ! ---------------------------- ! TODO: move these to a params structure? it=0 newstatbc=.true. ! This really shouldn't be here. ! General input per module ! ! This routine does need all processes, so not just xmaster ! Robert call all_input(par) ! Do check of params.txt to spot errors ! TODO: This shouldn't be in the main if (xmaster) call readkey('params.txt','checkparams',dummystring) ! TODO: We're not stepping into the timeloop just yet.... call writelog('ls','','Stepping into the time loop ....') ! writelog is xmaster aware #ifdef USEMPI call distribute_par(par) #endif if (xmaster) then call writelog('l','' ,'------------------------------------') call writelog('ls','','Building Grid and Bathymetry and....') call writelog('ls','','Distributing wave energy across the directional space ....') call writelog('l','', '------------------------------------') endif ! Grid and bathymetry ! ! grid_bathy will allocate x,y,xz,yz,xu,yv,xw,yw,zb,zb0 only ! on master process ! call space_alloc_scalars(sglobal) s => sglobal call grid_bathy(s,par) ! s%nx and s%ny are available now #ifdef USEMPI call xmpi_determine_processor_grid(s%nx,s%ny,par%mpiboundary,error) if(xmaster) then if (error==1) then call writelog('els','','Unknown mpi division ',par%mpiboundary) call halt_program else call writelog('ls','','processor grid: ',xmpi_m,' X ',xmpi_n) endif endif #endif ! initialize timesteps call timestep_init(par, tpar) if (xmaster) then ! Jump into subroutine readtide call readtide (s,par) !Ap 15/10 ! runs oonly on master wwvv call readwind (s,par) !Robert 8/7/2009 only on master call writelog('ls','','Initializing .....') ! Initialisations call flow_init (s,par) ! Always do this works only on master process call wave_init (s,par) ! Always do this wave_init only works on master process call gwinit(par,s) ! works only on master process call sed_init (s,par) ! works only on master process endif #ifdef USEMPI ! some of par has been changed, so: call distribute_par(par) #endif #ifdef USEMPI s => slocal ! ! determine how to divide the submatrices on the processor grid ! distribute all values in sglobal to slocal ! nx and ny will be adjusted in slocal ! arrays is,js,lm,ln (describing the distribution) will ! be filled in slocal ! Note: slocal is available on all nodes, including master ! call space_distribute_space(sglobal,slocal,par) !call space_consistency(slocal,'ALL') #endif ! Output means have to be initialized here after distribute space if (par%nmeanvar>0) call means_init(sglobal,slocal,par) call output_init(sglobal, slocal, par, tpar) ! update times at which we need output call outputtimes_update(par, tpar) ! Store first timestep (always) end function init integer(c_int) function outputext() bind(C, name="outputext") !DEC$ ATTRIBUTES DLLEXPORT::outputext outputext = 0 ! Calculate new output times, so we know when to stop call output(s,sglobal,par,tpar) end function outputext integer(c_int) function executestep() bind(C, name="executestep") !DEC$ ATTRIBUTES DLLEXPORT::executestep ! ---------------------------- ! This is the main time loop ! ---------------------------- executestep = 0 ! do while (par%t0) call flow_bc (s,par) !Dano moved here, after (long) wave bc generation if (it==0) then #ifdef USEMPI t01 = MPI_Wtime() #endif endif ! Wave timestep if (par%swave==1) then if (trim(par%instat) == 'stat' .or. trim(par%instat) == 'stat_table') then #ifdef USEMPI call wave_timestep(s,par) if ((abs(mod(par%t,par%wavint))<0.000001d0).or.newstatbc) then newstatbc=.false. endif #else if ((abs(mod(par%t,par%wavint))<0.000001d0).or.newstatbc) then call wave_stationary(s,par) newstatbc=.false. endif #endif else newstatbc=.false. call wave_timestep(s,par) endif endif ! Flow timestep if (par%gwflow==1) call gwflow(par,s) if (par%flow+par%nonh>0) call flow_timestep (s,par) if (par%ndrifter>0) call drifter (s,par) ! Suspended transport if(par%sedtrans==1) call transus(s,par) ! Bed level update if (par%morphology==1) call bed_update(s,par) ! Calculate new output times, so we know when to stop ! call output(s,sglobal,par,tpar) ! enddo end function executestep integer(c_int) function getdoubleparameter(name,value, length) bind(C,name="getdoubleparameter") !DEC$ ATTRIBUTES DLLEXPORT::getdoubleparameter USE iso_c_binding ! use inout otherwise things break real(c_double), intent(inout) :: value ! and we need the string length .... integer(c_int),value ,intent(in) :: length ! String character(kind=c_char),intent(in) :: name(length) ! Transform name to a fortran character... character(length) :: myname myname = char_array_to_string(name, length) select case (myname) case ('t') value = par%t case ('tstop') value = par%tstop case default value = -99.0d0 end select getdoubleparameter = 0 end function getdoubleparameter integer(c_int) function setdoubleparameter(name,value, length) bind(C,name="setdoubleparameter") !DEC$ ATTRIBUTES DLLEXPORT::setdoubleparameter USE iso_c_binding ! use inout otherwise things break real(c_double), intent(in) :: value ! and we need the string length .... integer(c_int), value ,intent(in) :: length ! String character(kind=c_char),intent(in) :: name(length) ! Transform name to a fortran character... character(length) :: myname myname = char_array_to_string(name, length) select case (myname) case ('t') par%t = value case ('tstop') par%tstop = value case ('tnext') tpar%tnext = value case default setdoubleparameter = -1 return end select setdoubleparameter = 0 end function setdoubleparameter integer(c_int) function getintparameter(name,value, length) bind(C,name="getintparameter") !DEC$ ATTRIBUTES DLLEXPORT::getintparameter USE iso_c_binding ! use inout otherwise things break integer(c_int), intent(inout) :: value ! and we need the string length .... integer(c_int),value ,intent(in) :: length ! String character(kind=c_char),intent(in) :: name(length) ! Transform name to a fortran character... character(length) :: myname myname = char_array_to_string(name, length) select case (myname) case ('nx') value = par%nx case ('ny') value = par%ny case default value = -99 end select getintparameter = 0 end function getintparameter integer(c_int) function getarray(name, x, length) bind(C, name="getarray") !DEC$ ATTRIBUTES DLLEXPORT::getarray ! use inout otherwise things break type(carraytype), intent(out) :: x ! and we need the string length .... integer(c_int),value ,intent(in) :: length ! String character(kind=c_char),intent(inout) :: name(length) character(length) :: myname integer :: index type(arraytype) :: array getarray = 1 myname = char_array_to_string(name, length) index = chartoindex(myname) call indextos(s,index,array) x = arrayf2c(array) getarray = 0 end function getarray integer(c_int) function get1ddoublearray(name, x, length) bind(C, name="get1ddoublearray") !DEC$ ATTRIBUTES DLLEXPORT::get1ddoublearray ! use inout otherwise things break type (c_ptr), value :: x ! and we need the string length .... integer(c_int),value ,intent(in) :: length ! String character(kind=c_char),intent(inout) :: name(length) character(length) :: myname integer :: index type(arraytype) :: array real(c_double), target, allocatable, dimension(:) :: r1 get1ddoublearray = -1 myname = char_array_to_string(name, length) index = chartoindex(myname) call indextos(s,index,array) allocate(r1(size(array%r1,1))) r1 = array%r1 ! array%r1 => r1 x = c_loc(r1) get1ddoublearray = 0 end function get1ddoublearray integer(c_int) function get2ddoublearray(name, x, length) bind(C, name="get2ddoublearray") !DEC$ ATTRIBUTES DLLEXPORT::get2ddoublearray ! use inout otherwise things break type (c_ptr), intent(out) :: x ! and we need the string length .... integer(c_int),value ,intent(in) :: length ! String character(kind=c_char),intent(in) :: name(length) character(length) :: myname integer :: index type(arraytype) :: array real(c_double), target, allocatable, dimension(:,:) :: r2 get2ddoublearray = -1 myname = char_array_to_string(name, length) index = chartoindex(myname) call indextos(s,index,array) allocate(r2(size(array%r2,1), size(array%r2,2))) r2 = array%r2 ! array%r2 => r2 x = c_loc(r2) get2ddoublearray = 0 end function get2ddoublearray integer(c_int) function finalize() bind(C, name="finalize") !DEC$ ATTRIBUTES DLLEXPORT::finalize ! ------------------------ ! End of main time loop ! ------------------------ ! ! Finish files finalize = -1 if(xmaster) then call cpu_time(tend) call writelog('ls','','Total calculation time: ',tend-tbegin,' seconds') #ifdef USEMPI t1 = MPI_Wtime() call writelog('ls','','Timing: procs: ',xmpi_size,' seconds: total:',t1-t0,'loop: ',t1-t01) #endif call writelog('ls','','End of program xbeach') endif call close_logfiles #ifdef USEMPI call xmpi_finalize #endif finalize = 0 end function finalize ! Utility functions function char_array_to_string(char_array, length) integer(c_int) :: length character(c_char) :: char_array(length) character(len=length) :: char_array_to_string integer :: i do i = 1, length char_array_to_string(i:i) = char_array(i) enddo end function end module libxbeach_module