236 investigate and separate current dependencies between observations module and the rest of the solver

.github/workflows/ubuntu.yml

@@ -29,6 +29,7 @@ jobs:
         mtln: ["ON", "OFF"]
         hdf: ["ON"]
         double-precision: ["OFF"]
+        new-output-module: ["OFF","ON"]
 
         include:
           # Disable by lack of space on github action
@@ -89,13 +90,15 @@ jobs:
           -DSEMBA_FDTD_ENABLE_MPI=${{matrix.mpi}} \
           -DSEMBA_FDTD_ENABLE_HDF=${{matrix.hdf}} \
           -DSEMBA_FDTD_ENABLE_MTLN=${{matrix.mtln}} \
-          -DSEMBA_FDTD_ENABLE_DOUBLE_PRECISION=${{matrix.double-precision}}
+          -DSEMBA_FDTD_ENABLE_DOUBLE_PRECISION=${{matrix.double-precision}} \
+          -DSEMBA_FDTD_ENABLE_OUTPUT_MODULE=${{matrix.new-output-module}}
         cmake --build build -j
 
     - name: Run unit tests 
       run: build/bin/fdtd_tests
 
     - name: Run python tests
+      if: matrix.new-output-module=='No'
       env:
         SEMBA_FDTD_ENABLE_MPI: ${{ matrix.mpi }}
         SEMBA_FDTD_ENABLE_MTLN: ${{ matrix.mtln }}

.gitmodules

@@ -20,4 +20,8 @@
 
 [submodule "external/googletest"]
 	path = external/googletest
-	url = https://github.com/google/googletest.git
+	url = https://github.com/google/googletest.git
+
+[submodule "external/VTKFortran"]
+	path = external/VTKFortran
+	url = https://github.com/szaghi/VTKFortran.git

CMakeLists.txt

@@ -16,7 +16,7 @@ option(SEMBA_FDTD_ENABLE_MPI "Use MPI" OFF)
 option(SEMBA_FDTD_ENABLE_HDF "Use HDF" ON)
 option(SEMBA_FDTD_ENABLE_MTLN "Use MTLN" ON)
 option(SEMBA_FDTD_ENABLE_SMBJSON "Use smbjson" ON)
-option(SEMBA_FDTD_ENABLE_DOUBLE_PRECISION "Use double precision (CompileWithReal8)" OFF)
+option(SEMBA_FDTD_ENABLE_DOUBLE_PRECISION "Use double precision (CompileWithReal8)" ON)
 
 option(SEMBA_FDTD_ENABLE_TEST "Compile tests" ON)
 
@@ -27,7 +27,21 @@ option(SEMBA_FDTD_EXECUTABLE "Compiles executable" ON)
 option(SEMBA_FDTD_MAIN_LIB "Compiles main library" ON)
 option(SEMBA_FDTD_COMPONENTS_LIB "Compiles components library" ON)
 option(SEMBA_FDTD_OUTPUTS_LIB "Compiles outputs library" ON)
+
+option(SEMBA_FDTD_ENABLE_OUTPUT_MODULE "Use new output module" OFF)
+
 # Compilation defines.
+if (CMAKE_SYSTEM_NAME STREQUAL "Windows")
+add_compile_definitions(_WIN32)
+elseif (CMAKE_SYSTEM_NAME STREQUAL "Darwin")
+add_compile_definitions(__APPLE__)
+elseif (CMAKE_SYSTEM_NAME STREQUAL "Linux")
+add_compile_definitions(__linux__)
+endif()
+
+if(SEMBA_FDTD_ENABLE_OUTPUT_MODULE)
+	add_definitions(-DCompileWithNewOutputModule)
+endif()
 if(SEMBA_FDTD_ENABLE_SMBJSON)
 	add_definitions(-DCompileWithSMBJSON)
 endif()
@@ -186,6 +200,13 @@ if (SEMBA_FDTD_ENABLE_MTLN)
 	endif()
 endif()
 
+add_subdirectory(src_utils)
+
+if (SEMBA_FDTD_ENABLE_OUTPUT_MODULE)
+	add_subdirectory(external/VTKFortran)
+	add_subdirectory(src_output)
+	set(OUTPUT_LIBRARIES fdtd-output)
+endif()
 add_subdirectory(src_conformal)
 set(CONFORMAL_LIBRARIES conformal)
 
@@ -194,6 +215,8 @@ if (SEMBA_FDTD_ENABLE_TEST)
 	add_subdirectory(test)
 endif()
 
+
+
 if(SEMBA_FDTD_COMPONENTS_LIB)
 	add_library(semba-components
 		"src_main_pub/anisotropic.F90"
@@ -256,7 +279,8 @@ if(SEMBA_FDTD_MAIN_LIB)
 		"src_main_pub/timestepping.F90"
 	)
 	target_link_libraries(semba-main
-		semba-outputs 
+		semba-outputs
+		${OUTPUT_LIBRARIES}
 		${SMBJSON_LIBRARIES} 
 		${MTLN_LIBRARIES})
 endif()

src_main_pub/fdetypes.F90

@@ -6,7 +6,6 @@
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 module  FDETYPES
 
-
 #ifdef CompileWithOpenMP
    use omp_lib
 #endif
@@ -186,6 +185,17 @@ module  FDETYPES
    integer (kind=4),  parameter  ::  iBloqueJx=100*iEx,iBloqueJy=100*iEy,iBloqueJz=100*iEz
    integer (kind=4),  parameter  ::  iBloqueMx=100*iHx,iBloqueMy=100*iHy,iBloqueMz=100*iHz
    !
+   integer (kind=4),  parameter :: VOLUMIC_M_MEASURE(3) = [iCur, iMEC, iMHC] !Module
+   integer (kind=4),  parameter :: VOLUMIC_X_MEASURE(3) = [iCurx, iExC, iHxC]
+   integer (kind=4),  parameter :: VOLUMIC_Y_MEASURE(3) = [iCury, iEyC, iHyC]
+   integer (kind=4),  parameter :: VOLUMIC_Z_MEASURE(3) = [iCurz, iEzC, iHzC]
+
+   integer (kind=4),  parameter :: MAGNETIC_FIELD_DIRECTION(3) = [iEx, iEy, iEz]
+   integer (kind=4),  parameter :: ELECTRIC_FIELD_DIRECTION(3) = [iHx, iHy, iHz]
+   integer (kind=4),  parameter :: CURRENT_MEASURE(4) = [iCur, iCurx, iCury, iCurz]
+   integer (kind=4),  parameter :: ELECTRIC_FIELD_MEASURE(4) = [iMEC, iExC, iEyC, iEzC]
+   integer (kind=4),  parameter :: MAGNETIC_FIELD_MEASURE(4) = [iMHC, iHxC, iHyC, iHzC]
+   !
    CHARACTER (LEN=*), PARAMETER  ::  SEPARADOR='______________'
    integer (kind=4), PARAMETER  ::  comi=1,fine=2, icoord=1,jcoord=2,kcoord=3
 
@@ -585,6 +595,7 @@ module  FDETYPES
 
 
    type  ::  MediaData_t
+      integer(kind=SINGLE) :: Id
       REAL (KIND=RKIND)          ::  Priority,Epr,Sigma,Mur,SigmaM
       logical :: sigmareasignado !solo afecta a un chequeo de errores en lumped 120123
       type (exists_t)            ::  Is
@@ -607,15 +618,15 @@ module  FDETYPES
       REAL (KIND=RKIND_tiempo)  ::  dt
       character (len=BUFSIZE) :: extraswitches
       !!
-      integer (kind=4)   ::  NumMedia,AllocMed
-      integer (kind=4)   ::  IniPMLMedia,EndPMLMedia
-      integer (kind=4)   ::  NumPlaneWaves,TimeSteps,InitialTimeStep
-      integer (kind=4)   ::  NumNodalSources
-      integer (kind=4)   ::  NumberRequest
+      integer (kind=SINGLE)   ::  NumMedia,AllocMed
+      integer (kind=SINGLE)   ::  IniPMLMedia,EndPMLMedia
+      integer (kind=SINGLE)   ::  NumPlaneWaves,TimeSteps,InitialTimeStep
+      integer (kind=SINGLE)   ::  NumNodalSources
+      integer (kind=SINGLE)   ::  NumberRequest = 0_SINGLE
       !!!
       REAL (KIND=RKIND)     , pointer, dimension ( : )        ::  LineX,LineY,LineZ
       REAL (KIND=RKIND)     , pointer, dimension ( : )        ::  DX,DY,DZ
-      integer (kind=4)                                        ::  AllocDxI,AllocDyI,AllocDzI,AllocDxE,AllocDyE,AllocDzE
+      integer (kind=SINGLE)                                        ::  AllocDxI,AllocDyI,AllocDzI,AllocDxE,AllocDyE,AllocDzE
       type (planeonde_t), pointer, dimension ( : )            ::  PlaneWave
       type (Border_t)                                         ::  Border
       type (PML_t)                                            ::  PML
@@ -632,6 +643,7 @@ module  FDETYPES
       logical  :: thereArePMLMagneticMedia
       CHARACTER (LEN=BUFSIZE) :: nEntradaRoot
       type (coorsxyzP)  ::  Punto
+
    end type
 
    type media_matrices_t
@@ -859,6 +871,8 @@ logical function direction_eq(a,b)
       direction_eq = direction_eq .and. (a%orientation == b%orientation)
 
    end function
+
+
 end module FDETYPES
 
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

src_main_pub/timestepping.F90

@@ -20,6 +20,10 @@ module Solver_mod
    use report
    use PostProcessing
    use Ilumina
+#ifdef CompileWithNewOutputModule
+   use output
+   use outputTypes
+#endif
    use Observa
    use BORDERS_other
    use Borders_CPML
@@ -1503,10 +1507,13 @@ subroutine initializeObservation()
          call MPI_Barrier(SUBCOMM_MPI,ierr)
 #endif
          write(dubuf,*) 'Init Observation...';  call print11(this%control%layoutnumber,dubuf)
+#ifdef CompileWithNewOutputModule
+         call init_outputs(this%sgg, this%media, this%sinPML_fullsize, this%bounds, this%control, this%thereAre%Observation, this%thereAre%wires)
+#else
          call InitObservation (this%sgg,this%media,this%tag_numbers, &
                                  this%thereAre%Observation,this%thereAre%wires,this%thereAre%FarFields,this%initialtimestep,this%lastexecutedtime, &
                                  this%sinPML_fullsize,this%eps0,this%mu0,this%bounds, this%control)
-
+#endif
          l_auxinput=this%thereAre%Observation.or.this%thereAre%FarFields
          l_auxoutput=l_auxinput
 
@@ -1772,6 +1779,10 @@ subroutine solver_run(this)
       real(kind=rkind), pointer, dimension (:,:,:) :: Ex, Ey, Ez, Hx, Hy, Hz
       real(kind=rkind), pointer, dimension (:) :: Idxe, Idye, Idze, Idxh, Idyh, Idzh, dxe, dye, dze, dxh, dyh, dzh
 
+#ifdef CompileWithNewOutputModule
+      type(fields_reference_t) :: fieldReference
+#endif
+
       logical :: call_timing, l_aux, flushFF, somethingdone, newsomethingdone
       integer :: i
       real (kind=rkind) :: pscale_alpha
@@ -1797,6 +1808,23 @@ subroutine solver_run(this)
 
       Idxe => this%Idxe; Idye => this%Idye; Idze => this%Idze; Idxh => this%Idxh; Idyh => this%Idyh; Idzh => this%Idzh; dxe => this%dxe; dye => this%dye; dze => this%dze; dxh => this%dxh; dyh => this%dyh; dzh => this%dzh
 
+#ifdef CompileWithNewOutputModule
+      fieldReference%E%x => this%Ex
+      fieldReference%E%y => this%Ey
+      fieldReference%E%z => this%Ez
+
+      fieldReference%E%deltax => this%dxe
+      fieldReference%E%deltay => this%dye
+      fieldReference%E%deltaz => this%dze
+
+      fieldReference%H%x => this%Hx
+      fieldReference%H%y => this%Hy
+      fieldReference%H%z => this%Hz
+
+      fieldReference%H%deltax => this%dxh
+      fieldReference%H%deltay => this%dyh
+      fieldReference%H%deltaz => this%dzh
+#endif
 
       ciclo_temporal :  DO while (this%n <= this%control%finaltimestep)
 
@@ -1872,9 +1900,11 @@ subroutine solver_run(this)
                       call print11(this%control%layoutnumber,SEPARADOR//separador//separador)
                       call print11(this%control%layoutnumber,dubuf)
                       call print11(this%control%layoutnumber,SEPARADOR//separador//separador)
-    !!
+#ifdef CompileWithNewOutputModule
+                      if (this%thereAre%Observation) call flush_outputs(this%sgg%tiempo, this%n, this%control, fieldReference, this%bounds, flushFF)
+#else
                       if (this%thereAre%Observation) call FlushObservationFiles(this%sgg,this%ini_save, this%n,this%control%layoutnumber, this%control%size, dxe, dye, dze, dxh, dyh, dzh,this%bounds,this%control%singlefilewrite,this%control%facesNF2FF,flushFF)
-                      !!
+#endif
 #ifdef CompileWithMPI
                       call MPI_Barrier(SUBCOMM_MPI,ierr)
 #endif
@@ -2023,11 +2053,21 @@ subroutine solver_run(this)
       subroutine updateAndFlush()
          integer(kind=4) :: mindum
          IF (this%thereAre%Observation) then
+#ifdef CompileWithNewOutputModule
+            if (this%n /= 0) then
+               call update_outputs(this%control, this%sgg%tiempo, this%n, fieldReference)
+               if (this%n>=this%ini_save+BuffObse)  then
+                  mindum=min(this%control%finaltimestep,this%ini_save+BuffObse)
+                  call flush_outputs(this%sgg%tiempo, this%n, this%control, fieldReference, this%bounds, .FALSE.)
+               endif
+            end if
+#else
             call UpdateObservation(this%sgg,this%media,this%tag_numbers, this%n,this%ini_save, Ex, Ey, Ez, Hx, Hy, Hz, dxe, dye, dze, dxh, dyh, dzh,this%control%wiresflavor,this%sinPML_fullsize,this%control%wirecrank, this%control%noconformalmapvtk,this%bounds)
             if (this%n>=this%ini_save+BuffObse)  then
                mindum=min(this%control%finaltimestep,this%ini_save+BuffObse)
                call FlushObservationFiles(this%sgg,this%ini_save,mindum,this%control%layoutnumber,this%control%size, dxe, dye, dze, dxh, dyh, dzh,this%bounds,this%control%singlefilewrite,this%control%facesNF2FF,.FALSE.) !no se flushean los farfields ahora
             endif
+#endif
          endif
       end subroutine
 
@@ -2726,12 +2766,35 @@ subroutine solver_end(this)
       logical :: dummylog, somethingdone, newsomethingdone
       character(len=bufsize) :: dubuf
 
+#ifdef CompileWithNewOutputModule
+      type(fields_reference_t) :: fieldReference
+#endif
+
+
 #ifdef CompileWithMPI
       integer (kind=4) :: ierr
 #endif
       Ex => this%Ex; Ey => this%Ey; Ez => this%Ez; Hx => this%Hx; Hy => this%Hy; Hz => this%Hz;
       dxe => this%dxe; dye => this%dye; dze => this%dze; dxh => this%dxh; dyh => this%dyh; dzh => this%dzh
 
+#ifdef CompileWithNewOutputModule
+      fieldReference%E%x => this%Ex
+      fieldReference%E%y => this%Ey
+      fieldReference%E%z => this%Ez
+
+      fieldReference%E%deltax => this%dxe
+      fieldReference%E%deltay => this%dye
+      fieldReference%E%deltaz => this%dze
+
+      fieldReference%H%x => this%Hx
+      fieldReference%H%y => this%Hy
+      fieldReference%H%z => this%Hz
+
+      fieldReference%H%deltax => this%dxh
+      fieldReference%H%deltay => this%dyh
+      fieldReference%H%deltaz => this%dzh
+#endif
+
 #ifdef CompileWithProfiling
       call nvtxEndRange
 #endif      
@@ -2777,8 +2840,12 @@ subroutine solver_end(this)
       call print11(this%control%layoutnumber,dubuf)
       call print11(this%control%layoutnumber,SEPARADOR//separador//separador)
       if (this%thereAre%Observation) THEN
+#ifdef CompileWithNewOutputModule
+         call flush_outputs(this%sgg%tiempo, this%n, this%control, fieldReference, this%bounds, .TRUE.)
+#else 
          call FlushObservationFiles(this%sgg,this%ini_save, this%n,this%control%layoutnumber, this%control%size, dxe, dye, dze, dxh, dyh, dzh,this%bounds,this%control%singlefilewrite,this%control%facesNF2FF,.TRUE.)
          call CloseObservationFiles(this%sgg,this%control%layoutnumber,this%control%size,this%control%singlefilewrite,this%initialtimestep,this%lastexecutedtime,this%control%resume) !dump the remaining to disk
+#endif
 #ifdef CompileWithMTLN      
          if (this%control%use_mtln_wires) then
             call FlushMTLNObservationFiles(this%control%nentradaroot, mtlnProblem = .false.)
@@ -2799,13 +2866,17 @@ subroutine solver_end(this)
       call MPI_Barrier(SUBCOMM_MPI,ierr)
 #endif
 
+#ifdef CompileWithNewOutputModule
+#else
+
       write(dubuf,'(a,i9)') 'INIT FINAL Postprocessing frequency domain probes, if any, at n= ',this%n
       call print11(this%control%layoutnumber,dubuf)
       write(dubuf,*) SEPARADOR//separador//separador
       call print11(this%control%layoutnumber,dubuf)
       somethingdone=.false.
       at=this%n*this%sgg%dt
       if (this%thereAre%Observation) call PostProcess(this%control%layoutnumber,this%control%size,this%sgg,this%control%nentradaroot,at,somethingdone,this%control%niapapostprocess,this%control%forceresampled)
+#endif
 #ifdef CompileWithMPI
       call MPI_Barrier(SUBCOMM_MPI,ierr)
       call MPI_AllReduce(somethingdone, newsomethingdone, 1_4, MPI_LOGICAL, MPI_LOR, SUBCOMM_MPI, ierr)
@@ -2831,7 +2902,6 @@ subroutine solver_end(this)
       somethingdone=.false.
 
       if (this%thereAre%Observation) call createvtk(this%control%layoutnumber,this%control%size,this%sgg,this%control%vtkindex,somethingdone,this%control%mpidir,this%media%sggMtag,this%control%dontwritevtk)
-
 #ifdef CompileWithMPI
       call MPI_Barrier(SUBCOMM_MPI,ierr)
       call MPI_AllReduce(somethingdone, newsomethingdone, 1_4, MPI_LOGICAL, MPI_LOR, SUBCOMM_MPI, ierr)

src_output/CMakeLists.txt

@@ -0,0 +1,19 @@
+add_library(fdtd-output 
+  "output.F90"
+  "outputTypes.F90"
+  "domain.F90"
+  "outputUtils.F90"
+  "volumicProbeUtils.F90"
+  "pointProbeOutput.F90"
+  "wireProbeOutput.F90"
+  "bulkProbeOutput.F90"
+  "movieProbeOutput.F90"
+  "frequencySliceProbeOutput.F90"
+  "farFieldProbeOutput.F90"
+)
+target_link_libraries(fdtd-output
+ semba-types
+ semba-components
+ fdtd-utils
+ VTKFortran::VTKFortran 
+)