Merge pull request #255 from OpenSEMBA/231-Unclear-porgram-flow-when-compiling-with-MTLN

231 Fix ambiguity in wire/multiwire in compilations with or without MTLN

Author: Alberto-o

doc/smbjson.md

@@ -51,7 +51,7 @@ The following entries are shared by several FDTD-JSON objects and have a common
 
 + `type` followed by a string, indicates the type of JSON object that. Some examples of types are `planewave` for `sources` objects, and `polyline` for `elements`.
 + `id` is a unique integer identifier for objects that belong to a list and which can be referenced by other objects. For instance, an element in the `elements` list must contain a `id` which can be referenced by a source in `sources` through its list of `elementIds`.
-+ `[name]` is an optional entry which is used to make the FDTD-JSON input human-readable, helping to identify inputs and outputs. Leading and trailing blank spaces are removed. Blank spaces are substituted by underscroes. The following characters are reserved and can't be used in a `name`: `@`.
++ `[name]` is an optional entry which is used to make the FDTD-JSON input human-readable, helping to identify inputs and outputs. Leading and trailing blank spaces are removed. Blank spaces are substituted by underscores. The following characters are reserved and can't be used in a `name`: `@`.
 
 ### `<general>`
 
@@ -60,9 +60,9 @@ This object must always be present and contains general information regarding th
 + `<timeStep>`: A real number indicating the time step used by the solver, in seconds. 
 + `<numberOfSteps>`: An integer for the number of steps which the solver will iterate.
 
-Addtionally, it may contain the following optional entry:
+Additionally, it may contain the following optional entry:
 
-+ `<mtlnProblem>` : A bool indicating whether the problem is a pure MTLN problem and will solved using only the MTLN solver. If it is not present, its default value is `false`
++ `<mtlnProblem>` : A boolean indicating whether the problem is a pure MTLN problem and will solved using only the MTLN solver. If it is not present, its default value is `false`
 + `<additionalArguments>` : A string with flags. Keep in mind that flags passed by console have higher priority. 
 
 **Example:**
@@ -78,7 +78,7 @@ Addtionally, it may contain the following optional entry:
 ### `[boundary]`
 This specifies the boundaries which will be used to terminate the computational domain. 
 If `boundary` is not present it defaults to a `mur` absorbing condition in all bounds.
-The entries within `boundary` are objects labelled with the place where they will be applied:
+The entries within `boundary` are objects labeled with the place where they will be applied:
 
 + `all`, or
 + `xLower`, `xUpper`, `yLower`, `yUpper`, `zLower` `zUpper`.
@@ -359,7 +359,7 @@ A `thinSlot` represents a gap between two conductive surfaces. Therefore it must
 
 ### `wire`
 
-A `wire`, or *thin wire*, represents an electrically conducting wire-like structure with a radius much smaller than the surrounding cell sizes. 
+A `wire`, or *thin wire*, represents an electrically conducting wire-like structure with a radius much smaller than the surrounding cell sizes. Materials of type `wire` can only be defined if the compilation flag for MTLN was OFF, i.e. `SEMBA_FDTD_ENABLE_MTLN = OFF`. Otherwise, the program execution will fail at runtime.
 These structures are solved by an algorithm similar to the one described in:
 
 ```
@@ -387,29 +387,14 @@ Materials of this type must contain:
 }
 ```
 
-A single wire might be surrounded by a dielectric material. In that case, the radius and the relative permittivity of the material are needed. 
-
-**Example:**
-
-```json
-{
-    "name": "WireWithDielectric",
-    "id": 2,
-    "type": "wire",
-    "radius": 0.0001,
-    "resistancePerMeter": 22.9e-3,
-    "dielectric" : {"radius": 0.001, "relativePermittivity" : 3}
-}
-```
-If the `dielectric` field is present but any of `radius` or `relativePermittivity` is absent, the parsing of the dielectric will fail.
 
 ### `shieldedMultiwire`
 
 A `shieldedMultiwire`, models $N+1$ electrical wires inside a bundled. The voltages and currents on these wires are solved by a multiconductor transmission lines (MTLN) solver described in:
 
     Paul, C. R. (2007). Analysis of multiconductor transmission lines. John Wiley & Sons.
 
-`shieldedMultiwire` materials are assumed to be contained within a `wire` or another `shieldedMultiwire` which is the external domain and is used as voltage reference. 
+`shieldedMultiwire` materials are assumed to be contained within an `unshieldedMultiwire` or another `shieldedMultiwire` which is the external domain and is used as voltage reference. Materials of type `shieldedMultiwire` and `unshieldedMultiwre` can only be defined if the compilation flag for MTLN was ON, i.e. `SEMBA_FDTD_ENABLE_MTLN = ON`. Otherwise, the program execution will fail at runtime.
 They must contain the following entries:
 
 + `<inductancePerMeter>` and `<capacitancePerMeter>` which must be matrices with a size $N \times N$. If the number of wires is equal to $1$, this property must be a $1 \times 1$ matrix, e.g `[[1e-7]]` 
@@ -449,7 +434,7 @@ They must contain the following entries:
 
 A `unshieldedMultiwire`, models a bundle of $N$ electrical wires. The charges and currents on these wires are solved using the model described in:
 
-    Berenger, J. P. A Multiwire formalism for the FDTD Method. IEEE Transactions on Electromagnetic Compatability. August, 2000.
+    Berenger, J. P. A Multiwire formalism for the FDTD Method. IEEE Transactions on Electromagnetic Compatibility. August, 2000.
 
 They must contain the following entries:
 
@@ -534,8 +519,7 @@ As with the rest of terminations, SPICE terminations have to be equivalents to 2
 
 ### `connector`
 
-The `connector` represents the physical connection of a bundle to a structure. `connector` assigns properties to the initial or last segment of a `wire`, a `shieldedMultiwire` or an `unshieldedMultiwire`. 
-This `wire` can be either a single wire or the outermost conductor of a `cable` bundle. The `connector` can have the following properties:
+The `connector` represents the physical connection of a bundle to a structure. `connector` assigns properties to the initial or last segment of a `wire`, a `shieldedMultiwire` or an `unshieldedMultiwire`. The `connector` can have the following properties:
 
 + `[resistances]`, an array of real numbers which will be converted to resistances per unit length and will replace the resistancePerMeter of that segment.
 + `[transferImpedancesPerMeter]`, an array of [transferImpedancePerMeter], as described in the [shieldedMultiwire](#shieldedMultiwire) section. 
@@ -689,7 +673,7 @@ In this example `elementId` points to a volume element, therefore `direction` mu
 }
 ```
 
-One important aspect to keep in mind when working with `bulkCurrent` with `electric field type` is its natural offset. This arises from the fact that to measure the electric current it is necessary to calculate the closed path integral of the magnetic field, then, the electric current is defined on the **dual mesh** of the inserted grid—i.e., the mesh corresponding to the magnetic field. The **dual mesh** is constructed by placing a point at the center of each cell in the original (primal) mesh and connecting these points.
+One important aspect to keep in mind when working with `bulkCurrent` with `electric field type` is its natural offset. This arises from the fact that to measure the electric current it is necessary to calculate the closed path integral of the magnetic field, then, the electric current is defined on the **dual mesh** of the inserted grid -- i.e., the mesh corresponding to the magnetic field. The **dual mesh** is constructed by placing a point at the center of each cell in the original (primal) mesh and connecting these points.
 
 Because of this, the code internally shifts the `bulkCurrent` you define to align with the dual mesh, causing a **half-cell offset**. In the case of surfaces, the coordinates **perpendicular** to the current flowing through the surface experience a **negative offset**, as shown in the figure below:
 
@@ -719,7 +703,7 @@ A `line` probe computes the electric field line integral along a given `polyline
 
 #### `farField`
 
-Probes of type `farField` perform a near to far field transformation of the electric and magnetic vector fields and are typically located in the scattered field region which is defined by a total/scatterd field excitation, e.g. [a planewave](#planewave). 
+Probes of type `farField` perform a near to far field transformation of the electric and magnetic vector fields and are typically located in the scattered field region which is defined by a total/scattered field excitation, e.g. [a planewave](#planewave). 
 They must be defined with a single `cell` element which must contain a single `interval` defining a cuboid. 
 The direction of  the radiated field $\hat{r}(\theta, \phi)$ is defined with `<theta>` and `<phi>`, which must contain `<initial>`, `<final>`, and `<step>`, expressed in degrees.
 The `domain` of a `farField` probe can only be of type `frequency`.

src_json_parser/smbjson.F90

@@ -183,12 +183,13 @@ function readProblemDescription(this) result (res)
       res%conformalRegs = this%readConformalRegions()
 
       ! Thin elements
+#ifdef CompileWithMTLN 
+      res%mtln = this%readMTLN()
+#else
       res%tWires = this%readThinWires()
+#endif
       res%tSlots = this%readThinSlots()
 
-#ifdef CompileWithMTLN
-      res%mtln = this%readMTLN()
-#endif
 
 
    end function
@@ -1868,10 +1869,17 @@ function buildBaseThinSlotComponent(cs) result(res)
    function readThinWires(this) result (res)
       class(parser_t) :: this
       type(ThinWires) :: res
-      type(materialAssociation_t), dimension(:), allocatable :: mAs
+      type(materialAssociation_t), dimension(:), allocatable :: mAs, mwires
       integer :: i, j
       logical :: found
 
+      mwires = this%getMaterialAssociations([ &
+                  J_MAT_TYPE_SHIELDED_MULTIWIRE//'  ',&
+                  J_MAT_TYPE_UNSHIELDED_MULTIWIRE    ])
+      if (size(mwires) /= 0) then 
+         call WarnErrReport('ERROR: shieldedMultiwires and unshieldedMultiwires can only be defined if compiled with MTLN', .true.)
+      end if
+
       mAs = this%getMaterialAssociations([J_MAT_TYPE_WIRE])
 
       ! Pre-allocates thin wires.
@@ -2510,21 +2518,25 @@ function readMTLN(this) result (mtln_res)
       class(parser_t) :: this
       type(mtln_t) :: mtln_res
       type(fhash_tbl_t) :: elemIdToPosition, elemIdToCable, connIdToConnector
-      type(materialAssociation_t), dimension(:), allocatable :: wires, multiwires, cables
+      type(materialAssociation_t), dimension(:), allocatable :: wires, cables
       class(cable_t), pointer :: ptr, read_cable
       integer :: i
 
+      wires = this%getMaterialAssociations([J_MAT_TYPE_WIRE])
+      if (size(wires) /=0 ) then
+            call WarnErrReport("ERROR: material type 'wire' is not allowed if compiled with MTLN", .true.)
+      end if
 
       cables = this%getMaterialAssociations([ &
                 J_MAT_TYPE_SHIELDED_MULTIWIRE//'  ',&
                 J_MAT_TYPE_UNSHIELDED_MULTIWIRE    ])
       ! spaces are needed to make strings have same length. 
       ! Why? Because of FORTRAN! It only accepts fixed length strings for arrays.
 
+
       mtln_res%connectors => readConnectors()
       call addConnIdToConnectorMap(connIdToConnector, mtln_res%connectors)
       if (size(cables) == 0) then 
-         mtln_res%has_multiwires = .false.
          mtln_res%time_step = 0
          mtln_res%number_of_steps = 0
          allocate(mtln_res%cables(0))
@@ -2533,7 +2545,6 @@ function readMTLN(this) result (mtln_res)
          return
       end if
 
-      mtln_res%has_multiwires = .true.
       mtln_res%time_step = this%getRealAt(this%root, J_GENERAL//'.'//J_GEN_TIME_STEP)
       mtln_res%number_of_steps = this%getRealAt(this%root, J_GENERAL//'.'//J_GEN_NUMBER_OF_STEPS)
 
@@ -2699,8 +2710,8 @@ function buildNetworks() result(res)
 
          allocate(aux_nodes(0))
          allocate(networks_coordinates(0))
-         cables = [ this%getMaterialAssociations([J_MAT_TYPE_WIRE]), &
-                    this%getMaterialAssociations([J_MAT_TYPE_UNSHIELDED_MULTIWIRE]), &
+         ! cables = [ this%getMaterialAssociations([J_MAT_TYPE_WIRE]), &
+         cables  = [this%getMaterialAssociations([J_MAT_TYPE_UNSHIELDED_MULTIWIRE]), &
                     this%getMaterialAssociations([J_MAT_TYPE_SHIELDED_MULTIWIRE]) ]
          do i = 1, size(cables)
             elemIds = cables(i)%elementIds

src_json_parser/smbjson_labels.F90

@@ -49,7 +49,6 @@ module smbjson_labels_mod
    character (len=*), parameter :: J_MAT_WIRE_REF_INDUCTANCE = "__referenceInductancePerMeter"
    character (len=*), parameter :: J_MAT_WIRE_DIELECTRIC = "dielectric"
    character (len=*), parameter :: J_MAT_WIRE_DIELECTRIC_RADIUS = "radius"
-   character (len=*), parameter :: J_MAT_WIRE_DIELECTRIC_PERMITTIVITY = "relativePermittivity"
 
    character (len=*), parameter :: J_MAT_LUMPED_MODEL = "model"
    character (len=*), parameter :: J_MAT_LUMPED_MODEL_RESISTOR = "resistor"

src_main_pub/fdetypes.F90

@@ -676,7 +676,7 @@ module  FDETYPES
                  fieldtotl,finishedwithsuccess, &
                  permitscaling,mtlnberenger,niapapostprocess, &
                  stochastic, verbose, dontwritevtk, &
-                 use_mtln_wires, resume_fromold, vtkindex,createh5bin,wirecrank,fatalerror
+                 resume_fromold, vtkindex,createh5bin,wirecrank,fatalerror
 #ifdef CompileWithConformal
       logical :: input_conformal_flag
 #endif

src_main_pub/interpreta_switches.F90

@@ -88,7 +88,6 @@ module interpreta_switches_m
          noconformalmapvtk, &
          createh5filefromsinglebin, &
          creditosyaprinteados, &
-         use_mtln_wires, &
          read_command_line
 
       integer(kind=4) :: &
@@ -696,8 +695,6 @@ subroutine interpreta(l, statuse)
                l%forceresampled = .true.
                l%opcionespararesumeo = trim(adjustl(l%opcionespararesumeo))//' '//trim(adjustl(l%chain))
 
-            case ('-mtlnwires')
-               l%use_mtln_wires = .true.
             CASE ('-wirethickness')
                i = i + 1
                CALL getcommandargument(l%chaininput, i, f, l%length, statuse, binaryPath)
@@ -1488,7 +1485,6 @@ subroutine print_help(l)
       CALL print11 (l%layoutnumber, '-wiresflavor {new/Slanted.or.experimental.or.slanted/transition/semistructured l%precision} : model for the wires    ')   
 #endif
       CALL print11(l%layoutnumber, '&                        (default '//trim(adjustl(l%wiresflavor))//')   ')
-      CALL print11(l%layoutnumber, '-mtlnwires             : Use mtln solver to advance wires currents ')
       CALL print11(l%layoutnumber, '-notaparrabos          : Do not remove extra double tails at the end of the wires ')
       CALL print11(l%layoutnumber, '&                        only available for the native format.             ')
       CALL print11(l%layoutnumber, '-intrawiresimplify     : Disable strict interpretation of .NFDE topology.  ')
@@ -1965,7 +1961,6 @@ subroutine default_flags(l)
       l%facesNF2FF%ab = .true.
       l%facesNF2FF%ar = .true.
       !defaults
-      l%use_mtln_wires = .false.
       l%read_command_line = .true.
       l%hay_slanted_wires = .false.
       l%forcing = .FALSE.

src_main_pub/observation.F90

@@ -4034,6 +4034,7 @@ subroutine FlushMTLNObservationFiles(nEntradaRoot, mtlnProblem)
 #endif
 
       mtln_solver => GetSolverPtr()
+      if (.not. allocated(mtln_solver%bundles)) return
       unit = 2000
       do i = 1, size(mtln_solver%bundles)
         do j = 1, size(mtln_solver%bundles(i)%probes)
@@ -4047,7 +4048,7 @@ subroutine FlushMTLNObservationFiles(nEntradaRoot, mtlnProblem)
             buffer = buffer//" "//"conductor_"//trim(adjustl(temp))
           end do
           write (unit, *) trim(buffer)
-          do k = 1, size(mtln_solver%bundles(i)%probes(j)%t)
+          do k = 1, mtln_solver%number_of_steps + 1
             buffer = ""
             write (temp, *) mtln_solver%bundles(i)%probes(j)%t(k)
             buffer = buffer//trim(temp)

src_main_pub/semba_fdtd.F90

@@ -646,6 +646,10 @@ subroutine semba_init(this, input_flags)
          CALL print11 (this%l%layoutnumber, '---> Conformal thin-gap solver: '//trim(adjustl(dubuf)))
          write(dubuf,*) this%l%run_with_dmma
          CALL print11 (this%l%layoutnumber, '---> DMMA thin-gap solver: '//trim(adjustl(dubuf)))
+#ifdef CompileWithMTLN
+         write(dubuf,'(a)') 'MTLN wires'
+         CALL print11 (this%l%layoutnumber, '---> Wire model: '//trim(adjustl(dubuf)))
+#else
          write(dubuf,'(a)') this%l%wiresflavor
          CALL print11 (this%l%layoutnumber, '---> Wire model: '//trim(adjustl(dubuf)))
          write(dubuf,'(a)') this%l%inductance_model
@@ -664,10 +668,11 @@ subroutine semba_init(this, input_flags)
          CALL print11 (this%l%layoutnumber, '---> Thin-wire double-tails removed: '//trim(adjustl(dubuf)))
          write(dubuf,*) this%l%fieldtotl                
          CALL print11 (this%l%layoutnumber, '---> Thin-wire -this%l%fieldtotl experimental switch: '//trim(adjustl(dubuf)))
+
          WRITE (dubuf,*) SEPARADOR // SEPARADOR // SEPARADOR
          CALL print11 (this%l%layoutnumber, dubuf)
+#endif
       endif
-
       IF (this%l%layoutnumber == 0) THEN
          call erasesignalingfiles(this%l%simu_devia)
       endif