XBeach
|
00001 module output_module 00002 00003 #ifdef HAVE_CONFIG_H 00004 #include "config.h" 00005 #endif 00006 use ncoutput_module, only: ncoutput, points_output_init, fortoutput_init 00007 #ifdef USENETCDF 00008 use ncoutput_module, only: ncoutput_init 00009 #endif 00010 use params, only: parameters 00011 use spaceparams 00012 use timestep_module, only: timepars, outputtimes_update 00013 use logging_module, only: writelog 00014 use xmpi_module 00015 use paramsconst 00016 ! IFDEF used in case netcdf support is not compiled, f.i. Windows (non-Cygwin) 00017 implicit none 00018 save 00019 private 00020 public output, output_error, output_init, log_progress 00021 00022 contains 00023 !_________________________________________________________________________________ 00024 00025 subroutine output_init(sglobal, slocal, par, tpar) 00026 implicit none 00027 type(spacepars), target, intent(inout) :: sglobal 00028 type(spacepars), target, intent(in) :: slocal 00029 type(parameters), intent(in) :: par 00030 type(timepars), intent(in) :: tpar 00031 00032 00033 ! get xz and yz in sglobal 00034 ! and initialize sglobal%collected and sglobal%precollected 00035 00036 #ifdef USEMPI 00037 sglobal%collected = .false. 00038 00039 call space_collect_mnem(sglobal,slocal,par,mnem_xz) 00040 00041 call space_collect_mnem(sglobal,slocal,par,mnem_yz) 00042 00043 if(xomaster) then 00044 sglobal%precollected = sglobal%collected 00045 endif 00046 #endif 00047 00048 ! initialize the correct output module (clean this up?, move to another module?) 00049 call points_output_init(sglobal,par) 00050 00051 select case(par%outputformat) 00052 00053 case(OUTPUTFORMAT_FORTRAN) 00054 ! only fortran 00055 call writelog('ls', '', 'Fortran outputformat') 00056 !call var_output_init(sglobal,slocal,par,tpar) 00057 call fortoutput_init(sglobal,par,tpar) 00058 case(OUTPUTFORMAT_NETCDF) 00059 ! only netcdf, stop if it's not build 00060 call writelog('ls', '', 'NetCDF outputformat') 00061 #ifdef USENETCDF 00062 call ncoutput_init(sglobal,slocal,par,tpar) 00063 #else 00064 call writelog('lse', '', 'This xbeach executable has no netcdf support. ', & 00065 'Rebuild with netcdf or run with outputformat=fortran') 00066 call halt_program 00067 #endif 00068 case(OUTPUTFORMAT_DEBUG) 00069 call writelog('ls', '', 'Debug outputformat, writing both netcdf and fortran output') 00070 #ifdef USENETCDF 00071 call writelog('ls', '', 'NetCDF outputformat') 00072 call ncoutput_init(sglobal,slocal,par,tpar) 00073 #endif 00074 !call var_output_init(sglobal,slocal,par,tpar) 00075 call writelog('ls', '', 'Fortran outputformat') 00076 call fortoutput_init(sglobal,par,tpar) 00077 endselect 00078 00079 end subroutine output_init 00080 !_________________________________________________________________________________ 00081 00082 subroutine output(sglobal,s,par,tpar,update) 00083 00084 use means_module 00085 use postprocessmod 00086 00087 implicit none 00088 00089 type(spacepars) :: s,sglobal 00090 type(parameters) :: par 00091 type(timepars) :: tpar 00092 00093 logical, optional :: update 00094 logical :: lupdate 00095 00096 logical :: end_program 00097 #ifdef USEMPI 00098 logical :: toall = .true. 00099 #endif 00100 00101 00102 if (present(update)) then 00103 lupdate = update 00104 else 00105 lupdate = .true. 00106 endif 00107 00108 ! update output times 00109 if(xcompute) then 00110 if (lupdate) call outputtimes_update(par, tpar) 00111 endif 00112 00113 ! update meanvars in current averaging period with current timestep 00114 if(par%nmeanvar .gt. 0) then 00115 if(xcompute) then 00116 if (par%t>tpar%tpm(1) .and. par%t<=tpar%tpm(size(tpar%tpm))) then 00117 call makeaverage(s,par) 00118 endif 00119 endif 00120 endif 00121 00122 ! write to log file is hotfile is generated 00123 if (tpar%outputh .and. par%writehotstart == 1) then 00124 call writelog('ls','','Writing hotstart file no. ',tpar%ith,' at t = ',par%t,' s') 00125 endif 00126 00127 end_program = .false. 00128 #ifdef USEMPI 00129 if(xcompute) then 00130 if(tpar%output) then 00131 call xmpi_send_sleep(xmpi_imaster,xmpi_omaster) ! wake up omaster 00132 call xmpi_bcast(end_program,toall) ! matching the xmpi_bcast 00133 ! ! in the do loop a few lines below 00134 call tell_xomaster_what_time_it_is ! matching the call a few lines below 00135 else 00136 return 00137 endif 00138 endif 00139 #endif 00140 00141 do 00142 ! xomaster will not leave this loop, except 00143 ! at the very end of the program 00144 #ifdef USEMPI 00145 if(xomaster) then 00146 call xmpi_send_sleep(xmpi_imaster,xmpi_omaster) 00147 call xmpi_bcast(end_program,toall) ! matching the xmpi_bcast 00148 ! ! above or the xmpi_bcast 00149 ! ! in finalize 00150 if(end_program) then 00151 call xmpi_barrier(toall) 00152 call xmpi_finalize 00153 stop 00154 endif 00155 call tell_xomaster_what_time_it_is ! matching the call a few lines above 00156 endif 00157 #endif 00158 00159 ! output 00160 00161 call ncoutput(sglobal,s,par, tpar) 00162 00163 ! clear averages after output of means 00164 if (tpar%outputm .and. tpar%itm>1) then 00165 call clearaverage(par) 00166 endif 00167 if (xcompute) exit 00168 enddo 00169 00170 #ifdef USEMPI 00171 contains 00172 subroutine tell_xomaster_what_time_it_is 00173 00174 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%tnext) 00175 00176 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%itg) 00177 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%itp) 00178 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%itm) 00179 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%itw) 00180 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%ith) 00181 00182 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%outputg) 00183 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%outputp) 00184 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%outputm) 00185 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%outputw) 00186 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%outputh) 00187 call xmpi_send(xmpi_imaster,xmpi_omaster,tpar%output) 00188 00189 call xmpi_send(xmpi_imaster,xmpi_omaster,par%t) 00190 00191 end subroutine tell_xomaster_what_time_it_is 00192 #endif 00193 end subroutine output 00194 !_________________________________________________________________________________ 00195 00196 subroutine output_error(s, sglobal, par, tpar) 00197 00198 use logging_module 00199 use xmpi_module, only: halt_program 00200 00201 implicit none 00202 00203 type(spacepars) :: s,sglobal 00204 type(parameters) :: par 00205 type(timepars) :: tpar 00206 00207 !call output(s, sglobal, par, tpar, update=.false.) 00208 call writelog('lse','','An extra output timestep is created to inquire the last timestep') 00209 call writelog('lse','',' before an error occured') 00210 call halt_program 00211 00212 end subroutine output_error 00213 00214 subroutine log_progress(par) 00215 00216 type(parameters) :: par 00217 logical,save :: firsttime = .true. 00218 integer,save :: day, ndt 00219 real*8,save :: tprev,percprev,sumdt,dtavg,t0 00220 real*8 :: tnow,percnow,tpredicted,tpredicted2,tpredmean 00221 integer,dimension(8) :: datetime 00222 00223 00224 if (firsttime) then 00225 day=0 00226 call date_and_time(VALUES=datetime) 00227 tprev = day*24.d0*3600.d0+datetime(5)*3600.d0+60.d0*datetime(6)+1.d0*datetime(7)+0.001d0*datetime(8) 00228 percprev = 0.d0 00229 t0 = tprev 00230 firsttime = .false. 00231 ndt = 0 00232 sumdt = 0.d0 00233 else 00234 if (par%timings .ne. 0) then 00235 ndt=ndt+1 00236 sumdt=sumdt+par%dt 00237 call date_and_time(VALUES=datetime) 00238 ! Current time in seconds 00239 tnow=day*24.d0*3600.d0+datetime(5)*3600.d0+60.d0*datetime(6)+1.d0*datetime(7)+0.001d0*datetime(8) 00240 ! Is it time to update the log again? 00241 if (tnow>=tprev+5.d0) then 00242 ! Percentage complete 00243 percnow = 100.d0*par%t/par%tstop 00244 ! Average time step in the last 5 seconds 00245 dtavg = sumdt/ndt 00246 call writelog('ls','(a,f5.1,a,f7.3,a)','Simulation ',percnow,' percent complete. Average dt ',sumdt/ndt,' seconds') 00247 ndt=0 00248 sumdt=0.d0 00249 ! Predict time based on percentage change rate in the last 5 seconds. 00250 tpredicted = 100.d0*(1.d0-par%t/par%tstop)/(max(percnow-percprev,0.01d0)/(tnow-tprev)) 00251 tpredicted2 = 100.d0*(1.d0-par%t/par%tstop)/(max(percnow,0.01d0)/(tnow-t0)) 00252 tpredmean = (tpredicted+tpredicted2)/2.d0 00253 ! Percentage complete: 00254 if (tpredmean>=3600) then 00255 call writelog('ls','(a,I3,a,I3,a)','Time remaining',& 00256 floor(tpredmean/3600.0d0),' hours and ',& 00257 nint((tpredmean-3600.0d0*floor(tpredmean/3600.0d0))/60.0d0),& 00258 ' minutes') 00259 elseif (tpredmean>=600) then 00260 call writelog('ls','(a,I3,a)','Time remaining ',& 00261 floor(tpredmean/60.0d0),' minutes') 00262 elseif (tpredmean>=60) then 00263 call writelog('ls','(a,I3,a,I3,a)','Time remaining ',& 00264 floor(tpredmean/60.0d0),' minutes and ',& 00265 nint((tpredmean-60.0d0*floor(tpredmean/60.0d0))),' seconds') 00266 else 00267 call writelog('ls','(a,I3,a)','Time remaining ',nint(tpredmean),' seconds') 00268 endif 00269 tprev=tnow 00270 percprev=percnow 00271 elseif (tnow<tprev-60.d0) then ! It's probably the next day 00272 day=day+1 00273 endif 00274 endif ! timings on 00275 endif ! firsttime logical 00276 end subroutine log_progress 00277 00278 end module output_module 00279