!! Copyright (C) Stichting Deltares, 2012-2017. !! !! This program is free software: you can redistribute it and/or modify !! it under the terms of the GNU General Public License version 3, !! as published by the Free Software Foundation. !! !! This program is distributed in the hope that it will be useful, !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! GNU General Public License for more details. !! !! You should have received a copy of the GNU General Public License !! along with this program. If not, see . !! !! contact: delft3d.support@deltares.nl !! Stichting Deltares !! P.O. Box 177 !! 2600 MH Delft, The Netherlands !! !! All indications and logos of, and references to registered trademarks !! of Stichting Deltares remain the property of Stichting Deltares. All !! rights reserved. !! Interface to access Delwaq from MT3D module MT3D logical :: MT3D_debug = .false. end module logical function MT3D_set_debug_status (set_MT3D_debug) !DEC$ ATTRIBUTES DLLEXPORT::MT3D_set_debug_status !DEC$ ATTRIBUTES DECORATE, ALIAS : 'MT3D_SET_DEBUG_STATUS' :: MT3D_set_debug_status use MT3D implicit none logical :: set_MT3D_debug MT3D_debug = set_MT3D_debug MT3D_set_debug_status = .true. end function logical function MT3D_processes_init (number_segments, number_substances, substances, & number_parameters, process_parameters, process_parameter_values, & number_processes, processes, deltt_in) !DEC$ ATTRIBUTES DLLEXPORT::MT3D_processes_init !DEC$ ATTRIBUTES DECORATE, ALIAS : 'MT3D_PROCESSES_INIT' :: MT3D_processes_init use waq_omi_priv use MT3D implicit none integer, intent(in) :: number_segments integer, intent(in) :: number_substances character(len=*), dimension(number_substances), intent(in) :: substances integer, intent(in) :: number_parameters character(len=*), dimension(number_parameters), intent(in) :: process_parameters real, dimension(number_segments, number_parameters), intent(in) :: process_parameter_values integer, intent(in) :: number_processes character(len=*), dimension(number_processes), intent(in) :: processes real, intent(in) :: deltt_in ! In days, real logical, external :: SetReferenceDate logical, external :: SetIntegrationOptions logical, external :: SetSimulationTimes logical, external :: DefineWQSchematisation logical, external :: DefineMonitoringLocations logical, external :: DefineWQDispersion logical, external :: DefineWQProcessesX logical, external :: DefineWQExtraOutputParameters logical, external :: SetCurrentValueScalarInit logical, external :: SetCurrentValueFieldInit logical, external :: SetOutputTimers logical, external :: SetInitialVolume integer, external :: ModelInitialize logical, external :: SetFlowData logical, external :: SetBoundaryConditions integer, parameter :: integration_method = 0 integer, dimension(1,4) :: pointers = reshape((/ 0, 0, 0, 0 /), (/ 1,4 /)) integer, dimension(4) :: number_exchanges = (/ 0, 0, 0, 0 /) real, dimension(3) :: dispersion = (/ 0.0, 0.0, 0.0 /) real, dimension(2) :: dispersion_length = (/ 1.0, 1.0 /) real, dimension(1) :: flow = (/ 0.0 /) real, dimension(1) :: area = (/ 1.0 /) real, allocatable :: substance_concentrations(:,:) real, allocatable :: substance_boundary_concentrations(:) real, allocatable :: volume(:) integer, allocatable :: substance_mult(:) real :: dvol_dt integer :: startTime = 0 !< Start time in seconds since the reference date/time integer :: endTime = huge(0) !< Stop time in seconds since the reference date/time integer :: deltt !< Time step in seconds character(len=8) :: cseg integer, dimension(number_segments) :: imon character(len=8), dimension(number_segments) :: cmon logical :: success, debug integer :: success_int, iseg, isub, ipar, i MT3D_processes_init = .false. deltt = int(deltt_in * 86400.0e0) ! From days (real) to seconds integer allocate (substance_boundary_concentrations(number_substances)) allocate (volume(number_segments)) allocate (substance_mult(number_substances)) substance_mult = 1 success = SetReferenceDate( 2014, 1, 1, 0, 0, 0 ) if(.not.success) return success = SetIntegrationOptions(integration_method, .true., .true., .true., .false., .false.) if(.not.success) return success = SetSimulationTimes( startTime, endTime, deltt ) if(.not.success) return success = DefineWQSchematisation( number_segments, pointers, number_exchanges ) if(.not.success) return if (mt3d_debug) then cseg = "s0000000" do iseg = 1, min(100,number_segments) imon(iseg) = iseg write(cseg(2:8), '(I7.7)') iseg cmon(iseg) = cseg end do success = DefineMonitoringLocations(imon, cmon, min(100,number_segments)) if(.not.success) return else success = DefineMonitoringLocations( (/0/), (/ ' ' /), 0 ) if(.not.success) return endif success = DefineWQDispersion(dispersion, dispersion_length ) ! dummy values if(.not.success) return success = DefineWQProcessesX( substances, substance_mult, number_substances, number_substances, & (/ 'DELT'/), 1 , & process_parameters, number_parameters, & processes, number_processes) if(.not.success) return success = DefineWQExtraOutputParameters( (/ 'DELT' /), 1 ) if(.not.success) return do isub = 1, number_substances success = SetCurrentValueScalarInit( substances(isub), 0.0 ) if(.not.success) return end do do ipar = 1, number_parameters success = SetCurrentValueFieldInit( process_parameters (ipar), process_parameter_values(1:number_segments,ipar)) if(.not.success) return end do if (mt3d_debug) then success = SetOutputTimers( 1, 0, endTime, deltt ) if(.not.success) return success = SetOutputTimers( 2, 0, endTime, deltt ) if(.not.success) return success = SetOutputTimers( 3, 0, endTime, deltt ) if(.not.success) return else success = SetOutputTimers( 1, endTime, endTime, endTime ) if(.not.success) return success = SetOutputTimers( 2, endTime, endTime, endTime ) if(.not.success) return success = SetOutputTimers( 3, endTime, endTime, endTime ) if(.not.success) return end if volume = 1.0e0 success = SetInitialVolume( volume ) if(.not.success) return success_int = ModelInitialize() if(.not.success_int.eq.0) return success = SetFlowData( volume, area, flow ) if(.not.success) return substance_boundary_concentrations = 0.0 success = SetBoundaryConditions( 1, substance_boundary_concentrations ) if(.not.success) return MT3D_processes_init = .true. end function MT3D_processes_init logical function MT3D_processes_setpar (number_parameters, number_segments, process_parameters, process_parameter_values) !DEC$ ATTRIBUTES DLLEXPORT::MT3D_processes_setpar !DEC$ ATTRIBUTES DECORATE, ALIAS : 'MT3D_PROCESSES_SETPAR' :: MT3D_processes_setpar ! use MT3D_PL implicit none integer, intent(in) :: number_parameters integer, intent(in) :: number_segments character(len=*), dimension(number_parameters), intent(in) :: process_parameters real, dimension(number_segments, number_parameters), intent(in) :: process_parameter_values logical, external :: SetCurrentValueFieldInit logical :: success integer :: ipar, i MT3D_processes_setpar = .false. do ipar = 1, number_parameters success = SetCurrentValueFieldInit( process_parameters (ipar), process_parameter_values(1:number_segments,ipar)) if(.not.success) return end do MT3D_processes_setpar = .true. end function MT3D_processes_setpar logical function MT3D_processes_step (number_segments, number_substances, substance_concentrations, deltt_in ) !DEC$ ATTRIBUTES DLLEXPORT::MT3D_processes_step !DEC$ ATTRIBUTES DECORATE, ALIAS : 'MT3D_PROCESSES_STEP' :: MT3D_processes_step ! use MT3D_PL use waq_omi_priv implicit none logical, external :: SetSimulationTimes integer, external :: ModelPerformTimeStep integer, intent(in) :: number_segments integer, intent(in) :: number_substances real, dimension(number_segments, number_substances), intent (inout) :: substance_concentrations real, intent(in) :: deltt_in ! In days integer :: startTime = 0 !< Start time in seconds since the reference date/time integer :: endTime !< Stop time in seconds since the reference date/time integer :: deltt !< Time step in seconds ! double precision, dimension(number_segments, number_substances) :: substance_concentrations_dbl logical :: success integer :: isub, i, success_int deltt = int(deltt_in * 86400.0e0) ! In seconds startTime = 0 endTime = huge(0) ! no end time until finalization success = SetSimulationTimes( startTime, endTime, deltt ) if(.not.success) return MT3D_processes_step = .false. do isub = 1, number_substances success = SetValuePriv( DLWQ_CONCENTRATION, isub, ODA_ALL_SEGMENTS, substance_concentrations(1:number_segments, isub), dlwq_set ) end do if(.not.success) return success_int = ModelPerformTimeStep() if(.not.success_int.eq.0) return do isub = 1, number_substances success = GetValuePriv( DLWQ_CONCENTRATION, isub, ODA_ALL_SEGMENTS, substance_concentrations(1:number_segments, isub) ) if(.not.success) return end do MT3D_processes_step = .true. end function MT3D_processes_step logical function MT3D_processes_fin () !DEC$ ATTRIBUTES DLLEXPORT::MT3D_processes_fin !DEC$ ATTRIBUTES DECORATE, ALIAS : 'MT3D_PROCESSES_FIN' :: MT3D_processes_fin implicit none integer, external :: ModelFinalize integer :: success_int success_int = ModelFinalize() MT3D_processes_fin = success_int.eq.0 end function MT3D_processes_fin