Rewriting surfaces and boundaries

src/modules/init/moduleInput.f90

@@ -48,6 +48,8 @@ MODULE moduleInput
       CALL readBoundaryEM(config)
       CALL checkStatus(config, "readBoundaryEM")
 
+      ! Read Physical Surfaces
+
       !Read Geometry
       CALL verboseError('Reading Geometry...')
       CALL readGeometry(config)
@@ -812,28 +814,120 @@ MODULE moduleInput
 
       END DO
 
-      ! Link physical surfaces to boundaries
-      call config%info('boundaries.particles.linking', found, n_children = nPhysicalSurfaces)
-      allocate(boundaryParticlesLinking(1:nPhysicalSurfaces))
-      allocate(speciesNames(1:nSpecies))
-      do b = 1, nPhysicalSurfaces
-        write(iString, '(i2)') b
-        object = 'boundary.particles.linking(' // trim(iString) // ')'
-        call config%get(object // '.physicalSurface', boundaryParticlesLinking(b)%physicalSurface, found) 
-        call config%get(object // '.model',           speciesNames,                                found) 
-        allocate(boundaryParticlesLinking(b)%speciesIndex(1:nSpecies))
-        do s = 1, nSpecies
-          boundaryParticlesLinking(b)%speciesIndex(s) = boundaryParticlesName_to_Index(speciesNames(s))
-
-        end do
-
-      end do
-
       !Init the list of particles from surfaces
       CALL OMP_INIT_LOCK(partSurfaces%lock)
 
     END SUBROUTINE readBoundaryParticle
 
+    SUBROUTINE readBoundaryEM(config)
+      USE moduleMesh
+      USE moduleOutput
+      USE moduleErrors
+      USE moduleEM
+      USE moduleSpecies
+      USE json_module
+      IMPLICIT NONE
+
+      TYPE(json_file), INTENT(inout):: config
+      CHARACTER(:), ALLOCATABLE:: object
+      LOGICAL:: found
+      INTEGER:: b
+      CHARACTER(2):: bString
+      character(len=100), allocatable:: modelName(:)
+
+      CALL config%info('boundaries.EM.models', found, n_children = nBoundariesEM)
+      
+      IF (found) THEN
+        ALLOCATE(boundariesEM(1:nBoundariesEM))
+
+      END IF
+      
+      do b = 1, nBoundaryEM
+        write(bString, '(I2)') b
+        object = 'boundaries.EM.models(' // TRIM(bString) // ')'
+
+        call boundariesEM(b)%init(config, object, b)
+
+      end do
+
+      ! TODO: Move this to the init of species
+      ALLOCATE(qSpecies(1:nSpecies))
+      DO s = 1, nSpecies
+        SELECT TYPE(sp => species(s)%obj)
+        TYPE IS (speciesCharged)
+          qSpecies(s) = sp%q
+
+        CLASS DEFAULT
+          qSpecies(s) = 0.D0
+
+        END SELECT
+
+      END DO
+
+    END SUBROUTINE readBoundaryEM
+
+    subroutine readPhysicalSurfaces(config)
+      use json_module
+      use moduleMesh
+      use moduleErrors
+      implicit none
+      
+      type(json_file), intent(inout):: config
+      character(:), allocatable:: object
+      logical:: found
+      integer:: ps
+      character(2):: psString
+      integer:: nParticleModels
+      character(len=100), allocatable:: particleModels(:)
+      character(len=100):: EMModel
+      integer:: s, boundaryIndex
+
+      call config%info('physicalSurfaces', found, n_children = nPhysicalSurfaces)
+      if (found) then
+        allocate(physicalSurfaces(1:nPhysicalSurfaces))
+
+      else
+        call criticalError('No physical surfaces found in the input file', 'readPhysicalSurfaces')
+
+      end if
+
+      do ps = 1, nPhysicalSurfaces
+        write(ps, '(I2)') ps
+        object = 'physicalSurfaces(' // TRIM(bString) // ')'
+
+        call config%get(object // '.index', physicalSurfaces(ps)%index, found) 
+        if (.not. found) then
+          call criticalError('Physical surface index not found', 'readPhysicalSurfaces')
+
+        end if
+        
+        ! Link models for particles
+        call config%info(object // '.particles', found, n_children = nParticleModels)
+        if ((.not. found) .or. &
+            (nParticleModels /= nSpecies)) then
+          call criticalError('Not enough models for particles provided', 'readPhysicalSurfaces')
+
+        end if
+        call config%get(object // '.particles', particleModels, found) 
+        allocate(physicalSurfaces(ps)%particles(1:nSpecies))
+        do s = 1, nSpecies
+          boundaryIndex = boundaryParticlesName_to_Index(particleModels(s))
+          physicalSurfaces(ps)%particles(s)%obj => boundaryParticles(boundaryIndex)%obj
+
+        end do
+
+        ! Link electromagnetic boundary condition
+        call config%get(object // '.EM', EMModel, found)
+        if (found) then
+          bouondaryIndex = boundaryEMName_to_Index(EMModel)
+          physicalSurfaces(ps)%EM => boundariesEM(boundaryIndex)%obj
+
+        end if
+
+      end do
+
+    end subroutine readPhysicalSurfaces
+
     !Read the geometry (mesh) for the case
     SUBROUTINE readGeometry(config)
       USE moduleMesh
@@ -1081,66 +1175,6 @@ MODULE moduleInput
 
     END SUBROUTINE readProbes
 
-    SUBROUTINE readBoundaryEM(config)
-      USE moduleMesh
-      USE moduleOutput
-      USE moduleErrors
-      USE moduleEM
-      USE moduleSpecies
-      USE json_module
-      IMPLICIT NONE
-
-      TYPE(json_file), INTENT(inout):: config
-      CHARACTER(:), ALLOCATABLE:: object
-      LOGICAL:: found
-      INTEGER:: b
-      CHARACTER(2):: bString
-      character(len=100), allocatable:: modelName(:)
-
-      CALL config%info('boundaries.EM.models', found, n_children = nBoundariesEM)
-      
-      IF (found) THEN
-        ALLOCATE(boundariesEM(1:nBoundariesEM))
-
-      END IF
-      
-      do b = 1, nBoundaryEM
-        write(bString, '(I2)') b
-        object = 'boundaries.EM.models(' // TRIM(bString) // ')'
-
-        call boundariesEM(b)%init(config, object, b)
-
-      end do
-
-      ! Link physical surfaces to boundaries
-      call config%info('boundaries.EM.linking', found, n_children = nPhysicalSurfaces)
-      allocate(boundaryEMLinking(1:nPhysicalSurfaces))
-      allocate(speciesNames(1:nSpecies))
-      do b = 1, nPhysicalSurfaces
-        write(iString, '(i2)') b
-        object = 'boundary.EM.linking(' // trim(iString) // ')'
-        call config%get(object // '.physicalSurface', boundaryEMLinking(b)%physicalSurface, found) 
-        call config%get(object // '.model',           modelName,                            found) 
-        boundaryEMLinking(b)%model => boundariesEM(boundaryEMName_to_Index(modelName(s)))%obj
-
-      end do
-
-      ! TODO: Move this to the init of species
-      ALLOCATE(qSpecies(1:nSpecies))
-      DO s = 1, nSpecies
-        SELECT TYPE(sp => species(s)%obj)
-        TYPE IS (speciesCharged)
-          qSpecies(s) = sp%q
-
-        CLASS DEFAULT
-          qSpecies(s) = 0.D0
-
-        END SELECT
-
-      END DO
-
-    END SUBROUTINE readBoundaryEM
-
     !Reads the injection of particles from the boundaries
     SUBROUTINE readInject(config)
       USE moduleSpecies

src/modules/mesh/makefile

@@ -23,6 +23,7 @@ moduleMesh.o: moduleMeshCommon.o moduleMesh.f90
 	$(FC) $(FCFLAGS) -c moduleMesh@elements.f90 -o $(OBJDIR)/$@
 	$(FC) $(FCFLAGS) -c moduleMesh@boundaryParticle.f90 -o $(OBJDIR)/$@
 	$(FC) $(FCFLAGS) -c moduleMesh@boundaryEM.f90 -o $(OBJDIR)/$@
+	$(FC) $(FCFLAGS) -c moduleMesh@surfaces.f90 -o $(OBJDIR)/$@
 
 inout.o: 3DCart.o 2DCyl.o 2DCart.o 1DRad.o 1DCart.o 0D.o
 	$(MAKE) -C inout all

src/modules/mesh/moduleMesh.f90

@@ -136,15 +136,14 @@ MODULE moduleMesh
 
   ABSTRACT INTERFACE
     !Inits the edge parameters
-    subroutine initEdge_interface(self, n, p, btPart, btEM)
+    subroutine initEdge_interface(self, n, p, physicalSurface)
       use moduleSpecies, only:nSpecies
       import:: meshEdge
 
       class(meshEdge), intent(out):: self
       integer, intent(in):: n
       integer, intent(in):: p(:)
-      integer, intent(in):: btPart(1:nSpecies)
+      integer, intent(in):: physicalSurface
-      integer, intent(in):: btEM
 
     end subroutine initEdge_interface
 
@@ -182,6 +181,12 @@ MODULE moduleMesh
 
   END TYPE meshEdgeCont
 
+  ! Array of pointers to edges.
+  type:: meshEdgePointer
+    class(meshEdge), pointer:: obj
+
+  end type meshEdgePointer
+
   !Parent of cell element
   TYPE, PUBLIC, ABSTRACT, EXTENDS(meshElement):: meshCell
     !Number of nodes in the cell
@@ -811,18 +816,10 @@ MODULE moduleMesh
 
   end type boundaryParticlePointer
 
-  type boundaryParticleLinking
-    integer:: physicalSurface
-    integer, allocatable, dimension(:):: speciesIndex
-
-  end type boundaryParticleLinking
-
   !Number of boundaries
-  INTEGER:: nBoundaries = 0, nPhysicalSurfaces = 0
+  INTEGER:: nBoundaries = 0
   !Array for boundaries
   TYPE(boundaryParticleCont), ALLOCATABLE, TARGET:: boundariesParticle(:)
-  !Array for linking boundaries
-  type(boundaryParticleLinking), allocatable, dimension(:):: boundariesParticleLinking
 
   ! BOUNDARY ELECTROMAGNETIC DEFINITIONS
   ! Generic type for electromagnetic boundary conditions
@@ -857,6 +854,12 @@ MODULE moduleMesh
 
     end subroutine addNodes
 
+    module function boundaryEMName_to_Index(boundaryName) result(bp)
+      character(:), allocatable:: boundaryName
+      integer:: bp
+
+    end function boundaryEMName_to_Index
+
   end interface
 
   abstract interface
@@ -879,6 +882,7 @@ MODULE moduleMesh
 
   end interface
 
+  ! Extended types
   TYPE, EXTENDS(boundaryEMGeneric):: boundaryEMDirichlet
     REAL(8):: potential
 
@@ -919,35 +923,34 @@ MODULE moduleMesh
 
   END TYPE boundaryEMCont
 
-  type:: boundaryEMLinking
-    integer:: physicalSurface
-    class(boundaryEMGeneric), pointer:: model => null()
-
-  end type boundaryEMLinking
-
   INTEGER:: nBoundariesEM
   TYPE(boundaryEMCont), ALLOCATABLE, target:: boundariesEM(:)
-  type(boundaryEMLinking), allocatable, dimension(:):: boundariesEMLinking
 
   !Information of charge and reference parameters for rho vector
   REAL(8), ALLOCATABLE:: qSpecies(:)
 
+  ! PHYSICAL SURFACES LINKING TO MESH ELEMENTS
+  ! Link physical surface to edges
+  type:: physicalSurface
+    integer:: index
+    class(meshNodePointer), allocatable:: nodes(:)
+    class(meshEdgePointer), allocatable:: edges(:)
+    class(boundaryParticlePointer), allocatable:: particles(:)
+    class(boundaryEMGeneric), pointer:: EM => null()
+
+  end type
+
+  integer:: nPhysicalSurfaces
+  type(physicalSurface), allocatable:: physicalSurfaces(:)
+
   ! Get ID from array based on physical surface
-  interface physicalSurface_to_id
+  interface
-    module function physicalSurface_to_boundaryParticleId(boundaryArray, physicalSurface) result(b)
+    module function physicalSurface_to_index(ps) result(index)
-      type(boundaryParticleLinking):: boundaryArray(:)
+      integer:: ps
-      integer:: physicalSurface
+      integer:: index
-      integer:: b
 
-    end function physicalSurface_to_boundaryParticleId
+    end function physicalSurface_to_index
 
-    module function physicalSurface_to_boundaryEMId(boundaryArray, physicalSurface) result(b)
+  end interface
-      type(boundaryEMLinking):: boundaryArray(:)
-      integer:: physicalSurface
-      integer:: b
-
-    end function physicalSurface_to_boundaryEMId
-
-  end interface physicalSurface_to_id
 
 END MODULE moduleMesh

src/modules/mesh/moduleMesh@boundaryEM.f90

@@ -214,25 +214,4 @@ submodule(moduleMesh) boundaryEM
 
     END SUBROUTINE applyDirichletTime
 
-    function physicalSurface_to_boundaryEMId(boundaryArray, physicalSurface) result(b)
-      implicit none
-
-      type(boundaryEMLinking):: boundaryArray(:)
-      integer:: physicalSurface
-      integer:: b
-      integer:: bt
-
-      b = 0
-      do bt=1, nPhysicalSurfaces
-        if (boundaryArray(bt)%physicalSurface == physicalSurface) then
-          b = boundaryArray(bt)%model%n
-          
-          exit
-
-        end if
-
-      end do
-
-    end function physicalSurface_to_boundaryEMId
-
 end submodule boundaryEM

src/modules/mesh/moduleMesh@boundaryParticle.f90

@@ -495,25 +495,4 @@ submodule(moduleMesh) boundaryParticle
 
     end subroutine genericTransparent
 
-    function physicalSurface_to_boundaryParticleId(boundaryArray, physicalSurface) result(b)
-      implicit none
-
-      type(boundaryParticleLinking):: boundaryArray(:)
-      integer:: physicalSurface
-      integer:: b
-      integer:: bt
-
-      b = 0
-      do bt=1, nPhysicalSurfaces
-        if (boundaryArray(bt)%physicalSurface == physicalSurface) then
-          b = bt
-          
-          exit
-
-        end if
-
-      end do
-
-    end function physicalSurface_to_boundaryParticleId
-
 end submodule boundaryParticle

src/modules/mesh/moduleMesh@surfaces.f90

@@ -0,0 +1,22 @@
+submodule(moduleMesh) surfaces
+  contains
+    function physicalSurface_to_index(ps) result(index)
+      implicit none
+
+      integer:: ps
+      integer:: index
+      integer:: i
+
+      index = 0
+      do i = 1, nPhysicalSurfaces
+        if (physicalSurfaces(i)%index == ps) then
+          index = i
+          exit
+
+        end if
+
+      end do
+
+    end function physicalSurface_to_index
+
+end submodule surfaces